1 DATA

#############################
#PERFORMANCE 
#############################

data_PERF <- import_data(dataset = "DATACOMPLET_PERF.csv", 
                         trait = "performance",
                         remove_testenvt = c("Grape","GF"), 
                         remove_pop = c("WT3"), 
                         remove_rate = NA)
## [1] "Data (541 and 602 tubes for the first and third generation respectively) where i) the number of eggs was NA (5 and 0 tubes for the first and third generation respectively); or ii) the number of adults  was NA (0 and 0 tubes for the first and third generation respectively); or iii) the number of eggs was zero -Emergence rate = NaN- (99 and 0 tubes for the first and third generation respectively); or iv) the number of adults was higher than the initial number of eggs (50 and 1 tubes for the first and third generation respectively) were not removed."
head(data_PERF)
##    Generation Experiment Original_environment   Population    Date_P Date_C_O
## 8          G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 13         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 14         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 16         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 20         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 24         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
##      Date_C_A Row Column Rack Test_environment Nb_eggs Obs_O Nb_adults Obs_A
## 8  19/10/2018  L1     C3    1           Cherry       3    LO         2    LO
## 13 19/10/2018  L1     C8    1       Blackberry       1    LO         2    LO
## 14 19/10/2018  L1     C9    1           Cherry       2    LO         1    LO
## 16 19/10/2018  L2     C1    1       Blackberry       1    LO         0    LO
## 20 19/10/2018  L2     C5    1       Strawberry       2    LO         2    LO
## 24 19/10/2018  L2     C9    1           Cherry       1    LO         1    LO
##    EggScore EggScoreFive EggScoreSmall SA IndicG0 IndicG2 SAIndicG0
## 8         1            1             1  0       1       0         0
## 13        1            1             1  1       1       0         1
## 14        1            1             1  0       1       0         0
## 16        1            1             1  1       1       0         1
## 20        1            1             1  0       1       0         0
## 24        1            1             1  0       1       0         0
##                fruit_hab             fruit_hab_ng     fruit_gen       hab_gen
## 8      Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
## 13 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0 Blackberry_G0
## 14     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
## 16 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0 Blackberry_G0
## 20 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0 Strawberry_G0
## 24     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
##            pop_gen      Rate
## 8  Blackberry45_G0 0.6666667
## 13 Blackberry45_G0 2.0000000
## 14 Blackberry45_G0 0.5000000
## 16 Blackberry45_G0 0.0000000
## 20 Blackberry45_G0 1.0000000
## 24 Blackberry45_G0 1.0000000
tapply(data_PERF$Nb_adults,list(data_PERF$Original_environment,data_PERF$Generation),length)
##             G0  G2
## Blackberry 169 320
## Cherry     233 143
## Strawberry  35 139
#Problem eggs

data_issue_eggs <- data_PERF[data_PERF$Nb_adults>data_PERF$Nb_eggs,]
#Ovservations with more eggs than edults
dim(data_issue_eggs)[1]
## [1] 51
#Ovservations with a difference equal or less than two eggs
sum(ifelse((data_issue_eggs$Nb_adults-data_issue_eggs$Nb_eggs)<=2,1,0))
## [1] 44
tapply(data_issue_eggs$Nb_eggs,data_issue_eggs$Test_environment,length)
## Blackberry     Cherry Strawberry 
##         18         17         16
tapply(data_issue_eggs$Nb_eggs,data_issue_eggs$Original_environment,length)
## Blackberry     Cherry Strawberry 
##         41          9          1
tapply(data_PERF$Nb_eggs,data_PERF$Original_environment,length)
## Blackberry     Cherry Strawberry 
##        489        376        174
mean(data_issue_eggs$Rate)
## [1] 2.44571
data_issue_eggs$Nb_eggs
##  [1]  1  1  2  1  1  1  1  1  1  3  3  1  2  1  5  1  1  3  2  1  2  3  2  1  2
## [26]  1  3  1  1  1  1  2  2  2  1  2  3  1  1  2  1  1  9  4  2 18  1  1  1  3
## [51] 78
data_issue_eggs$Nb_adults
##  [1]  2  3  4  3  3  2  2  2  2  5  5  2  4 12  7  2 13  5  3  2  3  5  6  2  3
## [26]  2  4  3  3  2  3  5  4  3  4  4  5  3  2  3  2  2 14  5  3 30  2  2  2  5
## [51] 80
data_issue_eggs$Rate
##  [1]  2.000000  3.000000  2.000000  3.000000  3.000000  2.000000  2.000000
##  [8]  2.000000  2.000000  1.666667  1.666667  2.000000  2.000000 12.000000
## [15]  1.400000  2.000000 13.000000  1.666667  1.500000  2.000000  1.500000
## [22]  1.666667  3.000000  2.000000  1.500000  2.000000  1.333333  3.000000
## [29]  3.000000  2.000000  3.000000  2.500000  2.000000  1.500000  4.000000
## [36]  2.000000  1.666667  3.000000  2.000000  1.500000  2.000000  2.000000
## [43]  1.555556  1.250000  1.500000  1.666667  2.000000  2.000000  2.000000
## [50]  1.666667  1.025641
#############################
#EMERGENCE RATE
#############################

#When Eggs>Adults: remove data
data_PERF_Rate_removed <- import_data(dataset = "DATACOMPLET_PERF.csv", 
                         trait = "performance",
                         remove_testenvt = c("Grape","GF"), 
                         remove_pop = c("WT3"), 
                         remove_rate = TRUE)
## [1] "Data (541 and 602 tubes for the first and third generation respectively) where i) the number of eggs was NA (5 and 0 tubes for the first and third generation respectively); or ii) the number of adults  was NA (0 and 0 tubes for the first and third generation respectively); or iii) the number of eggs was zero -Emergence rate = NaN- (99 and 0 tubes for the first and third generation respectively); or iv) the number of adults was higher than the number of eggs (50 and 1 tubes for the first and third generation respectively) were removed."
head(data_PERF_Rate_removed)
##    Generation Experiment Original_environment   Population    Date_P Date_C_O
## 8          G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 14         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 16         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 20         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 24         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
## 33         G0 Plasticity           Blackberry Blackberry45 3/10/2018  4/10/18
##      Date_C_A Row Column Rack Test_environment Nb_eggs Obs_O Nb_adults Obs_A
## 8  19/10/2018  L1     C3    1           Cherry       3    LO         2    LO
## 14 19/10/2018  L1     C9    1           Cherry       2    LO         1    LO
## 16 19/10/2018  L2     C1    1       Blackberry       1    LO         0    LO
## 20 19/10/2018  L2     C5    1       Strawberry       2    LO         2    LO
## 24 19/10/2018  L2     C9    1           Cherry       1    LO         1    LO
## 33 19/10/2018  L3     C8    1       Blackberry       2    LO         2    LO
##    EggScore EggScoreFive EggScoreSmall SA IndicG0 IndicG2 SAIndicG0
## 8         1            1             1  0       1       0         0
## 14        1            1             1  0       1       0         0
## 16        1            1             1  1       1       0         1
## 20        1            1             1  0       1       0         0
## 24        1            1             1  0       1       0         0
## 33        1            1             1  1       1       0         1
##                fruit_hab             fruit_hab_ng     fruit_gen       hab_gen
## 8      Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
## 14     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
## 16 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0 Blackberry_G0
## 20 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0 Strawberry_G0
## 24     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0     Cherry_G0
## 33 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0 Blackberry_G0
##            pop_gen      Rate
## 8  Blackberry45_G0 0.6666667
## 14 Blackberry45_G0 0.5000000
## 16 Blackberry45_G0 0.0000000
## 20 Blackberry45_G0 1.0000000
## 24 Blackberry45_G0 1.0000000
## 33 Blackberry45_G0 1.0000000
tapply(data_PERF_Rate_removed$Nb_adults,list(data_PERF_Rate_removed$Original_environment,
                                     data_PERF_Rate_removed$Generation),length)
##             G0  G2
## Blackberry 129 319
## Cherry     224 143
## Strawberry  34 139
tapply(data_PERF_Rate_removed$Rate,data_PERF_Rate_removed$Generation,mean)
##        G0        G2 
## 0.3593365 0.1857181
#When Eggs>Adults: Rate=1
data_PERF_Rate <- import_data(dataset = "DATACOMPLET_PERF.csv", 
                         trait = "performance",
                         remove_testenvt = c("Grape","GF"), 
                         remove_pop = c("WT3"), 
                         remove_rate = "Replace")
## [1] "Data (541 and 602 tubes for the first and third generation respectively) where i) the number of eggs was NA (5 and 0 tubes for the first and third generation respectively); or ii) the number of adults  was NA (0 and 0 tubes for the first and third generation respectively); or iii) the number of eggs was zero -Emergence rate = NaN- (99 and 0 tubes for the first and third generation respectively); or iv) the number of adults was higher than the initial number of eggs (50 and 1 tubes for the first and third generation respectively) were not removed but the emergence rate was REPLACED by 1."
tapply(data_PERF_Rate$Nb_adults,list(data_PERF_Rate$Original_environment,
                                     data_PERF_Rate$Generation),length)
##             G0  G2
## Blackberry 169 320
## Cherry     233 143
## Strawberry  35 139
tapply(data_PERF_Rate$Rate,data_PERF_Rate$Generation,mean)
##        G0        G2 
## 0.4326390 0.1870708
###########################
#PREFERENCE 
###########################
data_PREF <- import_data(dataset = "DATACOMPLET_PREF.csv", 
                         trait = "preference",
                         remove_testenvt = NA, 
                         remove_pop = c("WT3"), 
                         remove_rate = NA)

head(data_PREF)
##   Generation Experiment BoxID     Date_P Original_environment   Population Line
## 1         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 2         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 3         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 4         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 5         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    2
## 6         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    2
##   Column Test_environment Nb_eggs   Date_C_O Obs_O SA IndicG0 IndicG2 SAIndicG0
## 1      1        Cranberry       0 13/12/2018    CD  0       1       0         0
## 2      2              Fig       0 13/12/2018    CD  0       1       0         0
## 3      3        Raspberry       0 13/12/2018    CD  0       1       0         0
## 4      4         Rosehips       0 13/12/2018    CD  0       1       0         0
## 5      1             Kiwi       0 13/12/2018    CD  0       1       0         0
## 6      2       Strawberry       1 13/12/2018    CD  0       1       0         0
##               fruit_hab             fruit_hab_ng     fruit_gen       hab_gen
## 1  Blackberry_Cranberry  Blackberry_Cranberry_G0 Blackberry_G0  Cranberry_G0
## 2        Blackberry_Fig        Blackberry_Fig_G0 Blackberry_G0        Fig_G0
## 3  Blackberry_Raspberry  Blackberry_Raspberry_G0 Blackberry_G0  Raspberry_G0
## 4   Blackberry_Rosehips   Blackberry_Rosehips_G0 Blackberry_G0   Rosehips_G0
## 5       Blackberry_Kiwi       Blackberry_Kiwi_G0 Blackberry_G0       Kiwi_G0
## 6 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0 Strawberry_G0
##           pop_gen
## 1 Blackberry33_G0
## 2 Blackberry33_G0
## 3 Blackberry33_G0
## 4 Blackberry33_G0
## 5 Blackberry33_G0
## 6 Blackberry33_G0
tapply(data_PREF$Nb_eggs,list(data_PREF$Original_environment,data_PREF$Generation),length)
##              G0   G2
## Blackberry  696 1176
## Cherry     1200  624
## Strawberry  252  480
###########################
#PREFERENCE 3 fruits
###########################
levels_test<-levels(data_PREF$Test_environment)
levels_original<-levels(data_PREF$Original_environment)


data_PREF_three <- import_data(dataset = "DATACOMPLET_PREF.csv", 
                         trait = "preference",
                         remove_testenvt = usefun::outersect(levels_test, 
                                                             levels_original), 
                         remove_pop = c("WT3"), 
                         remove_rate = NA)

head(data_PREF_three)
##    Generation Experiment BoxID     Date_P Original_environment   Population
## 6          G0 Plasticity   860 19/09/2018           Blackberry Blackberry33
## 9          G0 Plasticity   860 19/09/2018           Blackberry Blackberry33
## 11         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33
## 15         G0 Plasticity   884 19/09/2018           Blackberry Blackberry33
## 16         G0 Plasticity   884 19/09/2018           Blackberry Blackberry33
## 20         G0 Plasticity   884 19/09/2018           Blackberry Blackberry33
##    Line Column Test_environment Nb_eggs   Date_C_O Obs_O SA IndicG0 IndicG2
## 6     2      2       Strawberry       1 13/12/2018    CD  0       1       0
## 9     3      1       Blackberry       0 13/12/2018    CD  1       1       0
## 11    3      3           Cherry       0 13/12/2018    CD  0       1       0
## 15    1      3           Cherry       0 13/12/2018    CD  0       1       0
## 16    1      4       Strawberry       0 13/12/2018    CD  0       1       0
## 20    2      4       Blackberry       0 13/12/2018    CD  1       1       0
##    SAIndicG0             fruit_hab             fruit_hab_ng     fruit_gen
## 6          0 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0
## 9          1 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0
## 11         0     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0
## 15         0     Blackberry_Cherry     Blackberry_Cherry_G0 Blackberry_G0
## 16         0 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0
## 20         1 Blackberry_Blackberry Blackberry_Blackberry_G0 Blackberry_G0
##          hab_gen         pop_gen
## 6  Strawberry_G0 Blackberry33_G0
## 9  Blackberry_G0 Blackberry33_G0
## 11     Cherry_G0 Blackberry33_G0
## 15     Cherry_G0 Blackberry33_G0
## 16 Strawberry_G0 Blackberry33_G0
## 20 Blackberry_G0 Blackberry33_G0
tapply(data_PREF_three$Nb_eggs,list(data_PREF_three$Original_environment,
                                    data_PREF_three$Generation),length)
##             G0  G2
## Blackberry 174 294
## Cherry     300 156
## Strawberry  63 120
resume_design<-tapply(as.factor(data_PREF_three$BoxID),list(data_PREF_three$Population,
                                  data_PREF_three$Generation),length)
mean(resume_design[,1], na.rm = TRUE)/3
## [1] 7.782609
mean(resume_design[,2], na.rm = TRUE)/3
## [1] 7.916667
resume_design<-tapply(data_PERF$Nb_eggs,list(data_PERF$Population,
                                  data_PERF$Generation),length)
mean(resume_design[,1], na.rm = TRUE)/3
## [1] 6.069444
mean(resume_design[,2], na.rm = TRUE)/3
## [1] 8.026667
dim(data_PREF_three)
## [1] 1107   21
dim(data_PREF)
## [1] 4428   21
head(data_PREF)
##   Generation Experiment BoxID     Date_P Original_environment   Population Line
## 1         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 2         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 3         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 4         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    1
## 5         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    2
## 6         G0 Plasticity   860 19/09/2018           Blackberry Blackberry33    2
##   Column Test_environment Nb_eggs   Date_C_O Obs_O SA IndicG0 IndicG2 SAIndicG0
## 1      1        Cranberry       0 13/12/2018    CD  0       1       0         0
## 2      2              Fig       0 13/12/2018    CD  0       1       0         0
## 3      3        Raspberry       0 13/12/2018    CD  0       1       0         0
## 4      4         Rosehips       0 13/12/2018    CD  0       1       0         0
## 5      1             Kiwi       0 13/12/2018    CD  0       1       0         0
## 6      2       Strawberry       1 13/12/2018    CD  0       1       0         0
##               fruit_hab             fruit_hab_ng     fruit_gen       hab_gen
## 1  Blackberry_Cranberry  Blackberry_Cranberry_G0 Blackberry_G0  Cranberry_G0
## 2        Blackberry_Fig        Blackberry_Fig_G0 Blackberry_G0        Fig_G0
## 3  Blackberry_Raspberry  Blackberry_Raspberry_G0 Blackberry_G0  Raspberry_G0
## 4   Blackberry_Rosehips   Blackberry_Rosehips_G0 Blackberry_G0   Rosehips_G0
## 5       Blackberry_Kiwi       Blackberry_Kiwi_G0 Blackberry_G0       Kiwi_G0
## 6 Blackberry_Strawberry Blackberry_Strawberry_G0 Blackberry_G0 Strawberry_G0
##           pop_gen
## 1 Blackberry33_G0
## 2 Blackberry33_G0
## 3 Blackberry33_G0
## 4 Blackberry33_G0
## 5 Blackberry33_G0
## 6 Blackberry33_G0
levels(data_PREF$Population)
##  [1] "Blackberry31" "Blackberry32" "Blackberry33" "Blackberry34" "Blackberry35"
##  [6] "Blackberry36" "Blackberry37" "Blackberry38" "Blackberry39" "Blackberry40"
## [11] "Blackberry43" "Blackberry44" "Blackberry45" "Cherry103"    "Cherry104"   
## [16] "Cherry3"      "Cherry47"     "Cherry50"     "Cherry51"     "Cherry52"    
## [21] "Cherry6"      "Cherry7"      "Strawberry42" "Strawberry44" "Strawberry53"

2 FECUNDITY

2.1 Exploration data

glmm::binomial.glmm()$family.glmm
## [1] "binomial.glmm"
tapply(data_PERF$Nb_eggs, list(data_PERF$Original_environment, data_PERF$Test_environment, data_PERF$Generation), mean, na.rm = TRUE)
## , , G0
## 
##            Blackberry   Cherry Strawberry
## Blackberry    6.52459 10.22034   7.306122
## Cherry       23.35443 19.77500  12.283784
## Strawberry   21.30769 33.60000  16.833333
## 
## , , G2
## 
##            Blackberry   Cherry Strawberry
## Blackberry   103.9720 121.8491   98.71963
## Cherry       131.7174 133.8400  103.51064
## Strawberry   119.1739 119.4894  111.56522
ggplot2::ggplot(data = data_PERF[data_PERF$Generation=="G0",], 
                aes(x = Test_environment, y = Nb_eggs, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PERF, 
                aes(x = Nb_eggs, fill = Original_environment)) +
  facet_wrap( ~ Generation) +
  geom_histogram(position="identity", alpha=0.5) +
  theme_LO_sober

2.2 Analysis

#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment,
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PERF)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[4,2])/(1/anova(m1)[4, 1]))
## [1] 1.971666
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[4, 1])) 
## [1] 0.2548912
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[5,2])/(1/anova(m2)[5, 1]))
## [1] 2.264412
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[5, 1]) )
## [1] 0.2294357
## F test for SA
(Fratio_NonGen <- (anova(m2)[4,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 1.369663
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[6, 1]) )
## [1] 0.3263818
## Compute R2 for SA 
## Compute R2 = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(rsqgen <- 1-anova(m2)[5, 3]/((anova(m2)[3, 2]+anova(m2)[5, 2])/(anova(m2)[3, 1]+anova(m2)[5, 1])))
## [1] 0.240181
(rsqng <- 1-anova(m2)[6, 3]/((anova(m2)[4, 2]+anova(m2)[6, 2])/(anova(m2)[4, 1]+anova(m2)[6, 1])))
## [1] 0.08459751
#######################################################
## Extract SA value   ###
#######################################################

#######################################################
## Extract SA value   ###
#######################################################
#Model used previously
m2 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF)
summary(m2)
##                                                Df Sum Sq Mean Sq F value
## pop_gen                                        48 1840.2   38.34  79.495
## hab_gen                                         4   15.6    3.90   8.085
## SA                                              1    2.7    2.73   5.655
## SA:IndicG0                                      1    1.0    0.99   2.054
## Original_environment:Test_environment           3    3.6    1.20   2.497
## IndicG0:Original_environment:Test_environment   3    2.2    0.72   1.499
## Residuals                                     978  471.7    0.48        
##                                                 Pr(>F)    
## pop_gen                                        < 2e-16 ***
## hab_gen                                       2.05e-06 ***
## SA                                              0.0176 *  
## SA:IndicG0                                      0.1522    
## Original_environment:Test_environment           0.0584 .  
## IndicG0:Original_environment:Test_environment   0.2132    
## Residuals                                                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Model with pophab interaction to have corret SA estimates
m3 <- lm(log(Nb_eggs+1)~ pop_gen + hab_gen + SA*IndicG0 + SA  , 
         data = data_PERF)

#Estimates SA
cf <-coef(summary(m3,complete = TRUE)) 
indic <- grep("SA",rownames(cf))
SAcoef <- cf[indic,1]
#names(SAcoef) <- c("SAGen", "SANonGen")
SAcoef
##         SA1 SA1:IndicG0 
##  -0.0575654  -0.1411427
#Estimates se
indic <- grep("SA",rownames(vcov(m3)))
seSAcoef <- sqrt(diag(vcov(m3)[indic,indic]))
#names(SAcoef) <- c("seSAGen", "seSANonGen")
seSAcoef
##         SA1 SA1:IndicG0 
##  0.06306995  0.09879136

##Plot

(PLOT_eggs_G0 <- plot_RTP_residuals(dataset = data_PERF, trait = "Nb_eggs", gen = "G0"))

(PLOT_eggs_G2 <- plot_RTP_residuals(dataset = data_PERF, trait = "Nb_eggs", gen = "G2"))

(PLOT_GEN_eggs_G0 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF, 
                                                   trait = "Nb_eggs", 
                                                   effect = "Non-genetic"))

(PLOT_GEN_eggs_G2 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF, 
                                                       trait = "Nb_eggs", 
                                                   effect = "Genetic"))

3 NUMBER OF ADULTS

3.1 Exploration data

tapply(data_PERF$Nb_adults, list(data_PERF$Original_environment, data_PERF$Test_environment, data_PERF$Generation), mean, na.rm = TRUE)
## , , G0
## 
##            Blackberry    Cherry Strawberry
## Blackberry   4.508197  4.457627   3.816327
## Cherry       6.873418  5.675000   2.378378
## Strawberry  12.230769 16.000000   5.416667
## 
## , , G2
## 
##            Blackberry   Cherry Strawberry
## Blackberry   27.31776 23.18868   13.15888
## Cherry       22.91304 23.32000   11.53191
## Strawberry   25.65217 22.36170   16.04348
ggplot2::ggplot(data = data_PERF[data_PERF$Generation=="G0",], 
                aes(x = Test_environment, y = Nb_adults, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PERF[data_PERF$Generation=="G2",], 
                aes(x = Test_environment, y = Nb_adults, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PERF, 
                aes(x = Nb_adults, fill = Original_environment)) +
  facet_wrap( ~ Generation) +
  geom_histogram(position="identity", alpha=0.5) +
  theme_LO_sober

3.2 Analysis

#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment,
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PERF)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[4,2])/(1/anova(m1)[4, 1]))
## [1] 2.01345
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[4, 1])) 
## [1] 0.2509626
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[5,2])/(1/anova(m2)[5, 1]))
## [1] 2.447544
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[5, 1]) )
## [1] 0.2156678
## F test for SA
(Fratio_NonGen <- (anova(m2)[4,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 1.111569
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[6, 1]) )
## [1] 0.3691494
## Compute R2 for SA 
## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(r2_SA_genet <- 1-(anova(m2)[5, 3]/((anova(m2)[3, 2]+anova(m2)[5, 2])/(anova(m2)[3, 1]+anova(m2)[5, 1]))))
## [1] 0.2657241
(r2_SA_nongenet <- 1-(anova(m2)[6, 3]/((anova(m2)[4, 2]+anova(m2)[6, 2])/(anova(m2)[4, 1]+anova(m2)[6, 1]))))
## [1] 0.02713541

3.3 Analysis with eggs as covariable

#######################################################
## Should we consider the number of eggs?   ###
#######################################################

# Original
m1 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF)


## Correction for number of eggs
m2 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            log(Nb_eggs+1), 
          data = data_PERF)


## With egg score
m3 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScore, 
          data = data_PERF)

## Compare with 5 egg scores
m4 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScoreFive, 
          data = data_PERF)

## Compare with EggScoreSmall
m5 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScoreSmall, 
          data = data_PERF)



MuMIn::model.sel(m1, m2, m3, m4, m5)
## Model selection table 
##       (Int) hab_gen pop_gen SA IG0:SA Org_env:Tst_env IG0:Org_env:Tst_env
## m2 -0.06648       +       +  +      +               +                   +
## m5  0.41840       +       +  +      +               +                   +
## m3  0.42840       +       +  +      +               +                   +
## m4  0.44140       +       +  +      +               +                   +
## m1  0.62560       +       +  +      +               +                   +
##    log(Nb_egg+1) EgS ESF ESS             family df    logLik   AICc  delta
## m2        0.5125             gaussian(identity) 63  -977.876 2090.0   0.00
## m5                         + gaussian(identity) 67 -1040.074 2223.5 133.51
## m3                 +         gaussian(identity) 65 -1045.327 2229.5 139.45
## m4                     +     gaussian(identity) 66 -1044.861 2230.8 140.80
## m1                           gaussian(identity) 62 -1118.184 2368.4 278.35
##    weight
## m2      1
## m5      0
## m3      0
## m4      0
## m1      0
## Models ranked by AICc(x)
# Models are not all fitted to the same data: because 6 tubes without Nb_eggs are missing for m2, m3, m4 and m5


## Cl= The best model is when the number of eggs is considered as a continuous variable
### model m2 provides a better description of the data than model m1 






#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment + log(Nb_eggs+1),
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PERF)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[5,2])/(1/anova(m1)[5, 1]))
## [1] 12.95037
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[5, 1])) 
## [1] 0.03679697
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(log(Nb_adults+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 + 
            log(Nb_eggs+1), 
          data = data_PERF)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 17.09331
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[6, 1]) )
## [1] 0.02567856
## F test for SA
(Fratio_NonGen <- (anova(m2)[5,2]/anova(m2)[7,2])/(1/anova(m2)[7, 1]))
## [1] 0.6257893
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[7, 1]) )
## [1] 0.4866763
## Compute R2 for SA 
## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(r2_SA_genet <- 1-(anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))))
## [1] 0.8009287
(r2_SA_nongenet <- 1-(anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))))
## [1] -0.1032081

##Plot

(PLOT_adult_G0 <- plot_RTP_residuals(dataset = data_PERF, trait = "Nb_adults", gen = "G0"))

(PLOT_adult_G2 <- plot_RTP_residuals(dataset = data_PERF, trait = "Nb_adults", gen = "G2"))

(PLOT_GEN_adult_G0 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF, 
                                                    trait = "Nb_adults", 
                                                    effect = "Non-genetic"))

(PLOT_GEN_adult_G2 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF, 
                                                    trait = "Nb_adults", 
                                                    effect = "Genetic"))

4 OFFPSRING PERFORMANCE

4.1 Exploration data

tapply(data_PERF_Rate$Rate, list(data_PERF_Rate$Original_environment, 
                                 data_PERF_Rate$Test_environment, 
                                 data_PERF_Rate$Generation), mean, na.rm = TRUE)
## , , G0
## 
##            Blackberry    Cherry Strawberry
## Blackberry  0.6951502 0.5960536  0.6508990
## Cherry      0.2708679 0.3454487  0.2493662
## Strawberry  0.3435592 0.3595403  0.3373855
## 
## , , G2
## 
##            Blackberry    Cherry Strawberry
## Blackberry  0.2648444 0.1923140  0.1345134
## Cherry      0.1855892 0.1801664  0.1103531
## Strawberry  0.2206513 0.1950778  0.1619433
ggplot2::ggplot(data = data_PERF_Rate[data_PERF_Rate$Generation=="G0",], 
                aes(x = Test_environment, y = Rate, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PERF_Rate[data_PERF_Rate$Generation=="G2",], 
                aes(x = Test_environment, y = Rate, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PERF_Rate, 
                aes(x = Rate, fill = Original_environment)) +
  facet_wrap( ~ Generation) +
  geom_histogram(position="identity", alpha=0.5) +
  theme_LO_sober

lattice::xyplot(Rate~Nb_eggs|Original_environment*Test_environment, 
                data=data_PERF_Rate)

lattice::xyplot(Rate~EggScore|Original_environment*Test_environment, 
                data=data_PERF_Rate)

lattice::bwplot(Rate~EggScore|Original_environment*Test_environment, 
                data=data_PERF_Rate)

lattice::bwplot(Rate~EggScoreFive|Original_environment*Test_environment, 
       data=data_PERF_Rate)

lattice::bwplot(Rate~EggScoreSmall|Original_environment*Test_environment, 
       data=data_PERF_Rate)

AvEmergenceRate <- tapply(data_PERF_Rate$Rate, 
                          list(data_PERF_Rate$EggScoreFive,
                              data_PERF_Rate$Original_environment,
                              data_PERF_Rate$Test_environment),mean)
tapply(data_PERF_Rate$Rate, list(data_PERF_Rate$EggScoreFive,
                                 data_PERF_Rate$Original_environment,
                                 data_PERF_Rate$Test_environment), length)
## , , Blackberry
## 
##   Blackberry Cherry Strawberry
## 1         71     72         14
## 2         43     22         15
## 3         41     17         20
## 4         10     10          8
## 5          3      4          2
## 
## , , Cherry
## 
##   Blackberry Cherry Strawberry
## 1         59     73          8
## 2         31     18         19
## 3         50     21         21
## 4         23     11          9
## 5          2      7         NA
## 
## , , Strawberry
## 
##   Blackberry Cherry Strawberry
## 1         55     75         18
## 2         55     24         13
## 3         39     16         16
## 4          7      6          9
## 5         NA     NA          2
AvEmergenceRate[, , "Blackberry"][, "Cherry"]/AvEmergenceRate[, , "Blackberry"][, "Blackberry"]
##         1         2         3         4         5 
## 0.4110203 0.9581939 0.7334304 0.9735309 0.1709919
AvEmergenceRate[, , "Cherry"][, "Blackberry"]/AvEmergenceRate[, , "Cherry"][, "Cherry"]
##         1         2         3         4         5 
## 1.6961620 0.8562518 1.0861981 0.8514873 1.4918117
AvEmergenceRate[, , "Blackberry"][, "Cherry"]/AvEmergenceRate[, , "Cherry"][, "Cherry"]
##         1         2         3         4         5 
## 0.7561587 1.0522488 1.0774054 0.9037047 0.4651315

4.2 Analysis without eggs

#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment,
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PERF_Rate)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[4,2])/(1/anova(m1)[4, 1]))
## [1] 20.07511
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[4, 1])) 
## [1] 0.02073057
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF_Rate)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[5,2])/(1/anova(m2)[5, 1]))
## [1] 30.92757
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[5, 1]) )
## [1] 0.01147033
## F test for SA
(Fratio_NonGen <- (anova(m2)[4,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 1.505315
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[6, 1]) )
## [1] 0.3073606
## Compute R2 for SA 
## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(r2_SA_genet <- 1-(anova(m2)[5, 3]/((anova(m2)[3, 2]+anova(m2)[5, 2])/(anova(m2)[3, 1]+anova(m2)[5, 1]))))
## [1] 0.8821018
(r2_SA_nongenet <- 1-(anova(m2)[6, 3]/((anova(m2)[4, 2]+anova(m2)[6, 2])/(anova(m2)[4, 1]+anova(m2)[6, 1]))))
## [1] 0.1121597

4.3 Variance of eggs as covariable

#Goal estimate: 
    # variance across tubes
    # variance across pop 
    # variance among host 

#### Mixte model
lm_vareggs <- lme4::lmer(log(Nb_eggs) ~ 1 +  (1|Population) + (1|Test_environment),
                         data = data_PERF_Rate)
summary(lm_vareggs)
## Linear mixed model fit by REML ['lmerMod']
## Formula: log(Nb_eggs) ~ 1 + (1 | Population) + (1 | Test_environment)
##    Data: data_PERF_Rate
## 
## REML criterion at convergence: 3855.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -2.8453 -0.5829  0.2357  0.7148  2.3291 
## 
## Random effects:
##  Groups           Name        Variance Std.Dev.
##  Population       (Intercept) 0.6633   0.8144  
##  Test_environment (Intercept) 0.0159   0.1261  
##  Residual                     2.2475   1.4992  
## Number of obs: 1039, groups:  Population, 25; Test_environment, 3
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)   3.3948     0.1856   18.29
(VAR_TESTENVT<-as.data.frame(lme4::VarCorr(lm_vareggs))$vcov[2])
## [1] 0.0159008
(VAR_POPULATION<-as.data.frame(lme4::VarCorr(lm_vareggs))$vcov[1])
## [1] 0.6633248
(VAR_TUBE<-as.data.frame(lme4::VarCorr(lm_vareggs))$vcov[3])
## [1] 2.247517
# Var host fruit: 0.01571
# Var population: 0.55243
# Var among tubes: 1.97954

## Extract variance components
VAR <- as.data.frame(lme4::VarCorr(lm_vareggs))$vcov

## Proportion of total variance due to variation among populations
(PropPop <- VAR[1]/sum(VAR))
## [1] 0.2266427
## Proportion of total variance  due to variation among fruits
(PropObs <- VAR[2]/sum(VAR))
## [1] 0.005432936
## Proportion of total variance due to within tube variance
(PropVarRes <- VAR[3]/sum(VAR)) 
## [1] 0.7679244
## Cl: High variance among tubes (78% of total variance)
## Cl: Substantial variation among popupation (22% of total variance)
## Cl: Low variance among hosts (>1% of total variance)


## ANOVA 
## Fit model with ANOVA
fit0 <- aov(log(Nb_eggs) ~ Population + Test_environment, data = data_PERF_Rate) 
anova(fit0)
## Analysis of Variance Table
## 
## Response: log(Nb_eggs)
##                    Df  Sum Sq Mean Sq F value  Pr(>F)    
## Population         24  611.54 25.4810 11.3416 < 2e-16 ***
## Test_environment    2   15.91  7.9528  3.5398 0.02938 *  
## Residuals        1012 2273.65  2.2467                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Variance among fruits 
(var_fruit <- (anova(fit0)[2,3]-anova(fit0)[3,3])/mean(tapply(data_PERF_Rate$Nb_eggs,
                                                              data_PERF_Rate$Test_environment,length)))
## [1] 0.01647582
#Variance among populations 
(var_pop <- (anova(fit0)[1,3]-anova(fit0)[3,3])/mean(tapply(data_PERF_Rate$Nb_eggs,
                                                            data_PERF_Rate$Population,length)))
## [1] 0.5590553
# Within group variance component = Residual variance
(var_tube <- (anova(fit0))[3,3])
## [1] 2.246687
sigma(fit0)^2 # the two match
## [1] 2.246687
#Plot 
mod1 <- lm(asin(sqrt(Rate)) ~ log(Nb_eggs) * Test_environment, data = data_PERF_Rate)

#Predict data
filldata <- data.frame(Nb_eggs = rep(seq(min(data_PERF_Rate$Nb_eggs),max(data_PERF_Rate$Nb_eggs)),3), 
                       Test_environment = rep(levels(data_PERF_Rate$Test_environment),
                                              each = length(seq(min(data_PERF_Rate$Nb_eggs),
                                                                max(data_PERF_Rate$Nb_eggs)))),
                       Estimates = NA)

filldata$Estimate_transformed <- predict(mod1, newdata=filldata, re.form=~0)
filldata$Estimates <- (sin(filldata$Estimate_transformed))^2


#If y=arcsin(sqrt(x)) then x=(sin(y))^2.

Plot_eggs <- ggplot2::ggplot(data = data_PERF_Rate, 
                aes(x = Nb_eggs, y = Rate, color = Test_environment)) +
  geom_point(size = 0.7)  +
  geom_line(data = filldata, aes(x = Nb_eggs, 
                                y = Estimates, 
                               colour = Test_environment), size = 0.9) + 
  scale_color_manual(name="Test fruit",   
                       breaks=c("Blackberry","Cherry","Strawberry"),
                       labels=c("Blackberry","Cherry","Strawberry"),
                       values=c("#301934","#BC3C6D", "#3FAA96"),
                       drop=FALSE) + 
  ylab("Offspring performance") +
  xlab("Number of eggs") +
  theme_LO_sober

 cowplot::save_plot(file =here::here("figures", "FigSX_DensityEggs_Effect.pdf"),
                   Plot_eggs,
                   base_height = 12/cm(1), base_width = 14/cm(1), dpi = 610)


 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 ######## Disucssion avec Nico 
var(log((data_PERF_Rate$Nb_adults+1)^2))
## [1] 6.260252
var(log((data_PERF_Rate$Nb_eggs+1)))
## [1] 2.251453
var(log((data_PERF_Rate$Nb_adults+1)))
## [1] 1.565063
var(log((data_PERF_Rate$Nb_eggs+1)))
## [1] 2.251453
var(log((data_PERF$Nb_adults)))
## [1] NaN
var(log((data_PERF$Nb_eggs)))
## [1] 2.794891
#The coefficient of variation (CV) is a statistical measure of the relative dispersion of data points in a data series around the mean. 

sd(data_PERF_Rate$Nb_adults[data_PERF_Rate$Generation=="G0"])/mean(data_PERF_Rate$Nb_adults[data_PERF_Rate$Generation=="G0"])
## [1] 1.768041
sd(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])/mean(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])
## [1] 1.376967
sd(data_PERF_Rate$Nb_adults[data_PERF_Rate$Generation=="G2"])/mean(data_PERF_Rate$Nb_adults[data_PERF_Rate$Generation=="G2"])
## [1] 0.821511
sd(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G2"])/mean(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G2"])
## [1] 0.385294
# 
# 
# sd(asin(sqrt(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G0"])))/mean(asin(sqrt(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G0"])))
# sd(asin(sqrt(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G2"])))/mean(asin(sqrt(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G2"])))
# 
# 
# sd(log(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"]))/mean(log(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"]))
# sd(log(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"]))/mean(log(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])) 



# 
# sd(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G0"])/mean(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G0"])
# sd(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G2"])/mean(data_PERF_Rate$Rate[data_PERF_Rate$Generation=="G2"])
# 
# 
# sd(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])/mean(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])
# sd(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])/mean(data_PERF_Rate$Nb_eggs[data_PERF_Rate$Generation=="G0"])


sd(data_PERF_Rate$Rate)/mean(data_PERF_Rate$Rate)
## [1] 1.01659
sd(data_PERF_Rate$Nb_eggs)/mean(data_PERF_Rate$Nb_eggs)
## [1] 0.8396281

4.4 Analysis with eggs as covariable

#######################################################
## Should we consider the number of eggs?   ###
#######################################################

# Original
m1 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0, 
          data = data_PERF_Rate)



## Correction for number of eggs
m2 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            log(Nb_eggs), 
          data = data_PERF_Rate)




## With egg score
m3 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScore, 
          data = data_PERF_Rate)

## Compare with 5 egg scores
m4 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScoreFive, 
          data = data_PERF_Rate)

## Compare with EggScoreSmall
m5 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            EggScoreSmall, 
          data = data_PERF_Rate)

## Correction for number of eggs
m6 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            Nb_eggs, 
          data = data_PERF_Rate)


data_PERF_Rate$Square_NbEggs <- data_PERF_Rate$Nb_eggs^2
## Correction for number of eggs
m7 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            Nb_eggs + Square_NbEggs, 
          data = data_PERF_Rate)


## Correction for number of eggs
m8 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            Square_NbEggs, 
          data = data_PERF_Rate)

MuMIn::model.sel(m1, m2, m3, m4, m5, m6, m7, m8)
## Model selection table 
##     (Int) hab_gen pop_gen SA IG0:SA Org_env:Tst_env IG0:Org_env:Tst_env
## m2 0.8300       +       +  +      +               +                   +
## m7 0.7880       +       +  +      +               +                   +
## m6 0.7427       +       +  +      +               +                   +
## m8 0.7305       +       +  +      +               +                   +
## m1 0.7400       +       +  +      +               +                   +
## m3 0.7471       +       +  +      +               +                   +
## m4 0.7460       +       +  +      +               +                   +
## m5 0.7474       +       +  +      +               +                   +
##    log(Nb_egg) EgS ESF ESS    Nb_egg    Sqr_NbE             family df   logLik
## m2     -0.1041                                  gaussian(identity) 63 -217.514
## m7                         -0.004514  1.424e-05 gaussian(identity) 64 -239.190
## m6                         -0.001382            gaussian(identity) 63 -246.327
## m8                                   -3.453e-06 gaussian(identity) 63 -251.902
## m1                                              gaussian(identity) 62 -254.885
## m3               +                              gaussian(identity) 65 -252.195
## m4                   +                          gaussian(identity) 66 -252.181
## m5                       +                      gaussian(identity) 67 -251.779
##     AICc delta weight
## m2 569.3  0.00      1
## m7 614.9 45.62      0
## m6 626.9 57.63      0
## m8 638.1 68.78      0
## m1 641.8 72.48      0
## m3 643.2 73.91      0
## m4 645.5 76.16      0
## m5 646.9 77.64      0
## Models ranked by AICc(x)
## Cl= The best model is when the number of eggs is considered as a continuous variable
### model m2 provides a better description of the data than model m1 







#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment +
            log(Nb_eggs),
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PERF_Rate)
  
  

## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[5,2])/(1/anova(m1)[5, 1]))
## [1] 16.94926
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[5, 1])) 
## [1] 0.02596669
## Compute R2 = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(rsqgen <- 1-anova(m1)[5, 3]/((anova(m1)[3, 2]+anova(m1)[4, 2])/(anova(m1)[3, 1]+anova(m1)[4, 1])))
## [1] 0.9842831
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            log(Nb_eggs), 
          data = data_PERF_Rate)


## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 18.92885
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[6, 1]) )
## [1] 0.02242804
## F test for SA
(Fratio_NonGen <- (anova(m2)[5,2]/anova(m2)[7,2])/(1/anova(m2)[7, 1]))
## [1] 0.7229752
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[7, 1]) )
## [1] 0.4576467
# Compute R2 = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
rsqgen <- 1-anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))
rsqng <- 1-anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))
  

## Compute R2 for SA 
## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(r2_SA_genet <- 1-(anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))))
## [1] 0.8175919
(r2_SA_nongenet <- 1-(anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))))
## [1] -0.07440953
#######################################################
## Extract SA value   ###
#######################################################
#Model used previously
m2 <- aov(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            log(Nb_eggs), 
          data = data_PERF_Rate)
summary(m2)
##                                                Df Sum Sq Mean Sq F value
## pop_gen                                        48  83.07   1.731  18.286
## hab_gen                                         4   2.24   0.561   5.925
## SA                                              1   1.06   1.062  11.222
## log(Nb_eggs)                                    1   6.91   6.912  73.030
## SA:IndicG0                                      1   0.08   0.078   0.820
## Original_environment:Test_environment           3   0.17   0.056   0.593
## IndicG0:Original_environment:Test_environment   3   0.32   0.107   1.134
## Residuals                                     977  92.47   0.095        
##                                                 Pr(>F)    
## pop_gen                                        < 2e-16 ***
## hab_gen                                       0.000103 ***
## SA                                            0.000839 ***
## log(Nb_eggs)                                   < 2e-16 ***
## SA:IndicG0                                    0.365390    
## Original_environment:Test_environment         0.619771    
## IndicG0:Original_environment:Test_environment 0.334151    
## Residuals                                                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#Model with pophab interaction to have corret SA estimates
m3 <- lm(asin(sqrt(Rate)) ~ pop_gen + hab_gen + SA*IndicG0 + SA  + log(Nb_eggs), 
         data = data_PERF_Rate)

#Estimates SA
cf <-coef(summary(m3,complete = TRUE)) 
indic <- grep("SA",rownames(cf))
SAcoef <- cf[indic,1]
#names(SAcoef) <- c("SAGen", "SANonGen")
SAcoef
##         SA1 SA1:IndicG0 
##  0.04226139  0.03955350
#Estimates se
indic <- grep("SA",rownames(vcov(m3)))
seSAcoef <- sqrt(diag(vcov(m3)[indic,indic]))
#names(SAcoef) <- c("seSAGen", "seSANonGen")
seSAcoef
##         SA1 SA1:IndicG0 
##  0.02785344  0.04366013

4.5 Heterogeneity across populations

data_Rate_G2 <- data_PERF_Rate[data_PERF_Rate$Generation=="G2",]
data_Rate_G2 <- data_Rate_G2[complete.cases(data_Rate_G2$Rate), ]

  
m0 <- lm(Rate~Original_environment*Test_environment, data=data_Rate_G2,)

summary(m0)
## 
## Call:
## lm(formula = Rate ~ Original_environment * Test_environment, 
##     data = data_Rate_G2)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.26484 -0.09978 -0.02467  0.07237  0.73516 
## 
## Coefficients:
##                                                           Estimate Std. Error
## (Intercept)                                                0.26484    0.01330
## Original_environmentCherry                                -0.07926    0.02426
## Original_environmentStrawberry                            -0.04419    0.02426
## Test_environmentCherry                                    -0.07253    0.01886
## Test_environmentStrawberry                                -0.13033    0.01881
## Original_environmentCherry:Test_environmentCherry          0.06711    0.03385
## Original_environmentStrawberry:Test_environmentCherry      0.04696    0.03421
## Original_environmentCherry:Test_environmentStrawberry      0.05509    0.03418
## Original_environmentStrawberry:Test_environmentStrawberry  0.07162    0.03431
##                                                           t value Pr(>|t|)    
## (Intercept)                                                19.909  < 2e-16 ***
## Original_environmentCherry                                 -3.267 0.001150 ** 
## Original_environmentStrawberry                             -1.822 0.069020 .  
## Test_environmentCherry                                     -3.846 0.000133 ***
## Test_environmentStrawberry                                 -6.928 1.12e-11 ***
## Original_environmentCherry:Test_environmentCherry           1.982 0.047891 *  
## Original_environmentStrawberry:Test_environmentCherry       1.373 0.170346    
## Original_environmentCherry:Test_environmentStrawberry       1.612 0.107534    
## Original_environmentStrawberry:Test_environmentStrawberry   2.088 0.037265 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1376 on 593 degrees of freedom
## Multiple R-squared:  0.1044, Adjusted R-squared:  0.09231 
## F-statistic:  8.64 on 8 and 593 DF,  p-value: 3.499e-11
tapply(data_Rate_G2$Rate, 
       list(data_Rate_G2$Original_environment, 
            data_Rate_G2$Test_environment), mean)
##            Blackberry    Cherry Strawberry
## Blackberry  0.2648444 0.1923140  0.1345134
## Cherry      0.1855892 0.1801664  0.1103531
## Strawberry  0.2206513 0.1950778  0.1619433
tapply(data_Rate_G2$Rate, 
       list(data_Rate_G2$Original_environment, 
            data_Rate_G2$Test_environment), var)
##            Blackberry      Cherry Strawberry
## Blackberry 0.03491514 0.013749525 0.01436498
## Cherry     0.02922031 0.019467766 0.01399160
## Strawberry 0.01555929 0.009464032 0.01139614
range(data_Rate_G2$Nb_adults, na.rm = TRUE)
## [1]   0 108
## Check number of eggs and adults
tapply(data_Rate_G2$Nb_eggs, 
       list(data_Rate_G2$Original_environment, 
            data_Rate_G2$Test_environment), mean)
##            Blackberry   Cherry Strawberry
## Blackberry   103.9720 121.8491   98.71963
## Cherry       131.7174 133.8400  103.51064
## Strawberry   119.1739 119.4894  111.56522
tapply(data_Rate_G2$Nb_adults, 
       list(data_Rate_G2$Original_environment, 
            data_Rate_G2$Test_environment), mean)
##            Blackberry   Cherry Strawberry
## Blackberry   27.31776 23.18868   13.15888
## Cherry       22.91304 23.32000   11.53191
## Strawberry   25.65217 22.36170   16.04348
## Check for the presence of negative correlations
m1 <- lm(asin(sqrt(Rate)) ~ Population + Test_environment,
         data=data_Rate_G2)
data_Rate_G2$resid <- residuals(m1)

meanbypopbytestenv <- as.data.frame(tapply(data_Rate_G2$resid, 
                                           list(data_Rate_G2$Population,
                                                data_Rate_G2$Test_environment), mean))


## Cherry ~ Blackberry
plot(meanbypopbytestenv$Blackberry, 
     meanbypopbytestenv$Cherry, 
     col=as.numeric(data_Rate_G2$Original_environment[match(rownames(meanbypopbytestenv),
                                                    data_Rate_G2$Population)]), 
     xlab="Blackberry", ylab="Cherry", pch=16)
legend("topright", levels(data_Rate_G2$Original_environment), 
       col=as.numeric(as.factor(levels(data_Rate_G2$Original_environment))), pch=16)

## Strawberry ~ Cherry
plot(meanbypopbytestenv$Cherry,  meanbypopbytestenv$Strawberry, col=as.numeric(data_Rate_G2$Original_environment[match(rownames(meanbypopbytestenv), data_Rate_G2$Population)]), xlab="Cherry", ylab="Strawberry", pch=16)
legend("bottomright", levels(data_Rate_G2$Original_environment), col=as.numeric(as.factor(levels(data_Rate_G2$Original_environment))), pch=16)

## Strawberry ~ Blackberry
plot(meanbypopbytestenv$Blackberry, meanbypopbytestenv$Strawberry, col=as.numeric(data_Rate_G2$Original_environment[match(rownames(meanbypopbytestenv), data_Rate_G2$Population)]), xlab="Blackberry", ylab="Strawberry", pch=16)
legend("topright", levels(data_Rate_G2$Original_environment), col=as.numeric(as.factor(levels(data_Rate_G2$Original_environment))), pch=16)

##Plot

(PLOT_rate_G0 <- plot_RTP_residuals(dataset = data_PERF_Rate, trait = "Rate", gen = "G0"))

(PLOT_rate_G2 <- plot_RTP_residuals(dataset = data_PERF_Rate, trait = "Rate", gen = "G2"))

(PLOT_GEN_rate_G0 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF_Rate, 
                                                    trait = "Rate", 
                                                    effect = "Non-genetic"))

(PLOT_GEN_rate_G2 <- plot_Genetic_Nongenetic_residuals(dataset = data_PERF_Rate, 
                                                    trait = "Rate", 
                                                    effect = "Genetic"))

5 OVIOSITION PREFERENCE

5.1 Twelve fruits

5.1.1 Exploration data

tapply(data_PREF$Nb_eggs, list(data_PREF$Original_environment, 
                                 data_PREF$Test_environment, 
                                 data_PREF$Generation), mean, na.rm = TRUE)
## , , G0
## 
##              Apricot Blackberry Blackcurrant    Cherry Cranberry       Fig
## Blackberry 0.2758621  0.7931034    0.1379310 0.4137931 0.2068966 0.2586207
## Cherry     0.3300000  0.5000000    1.2600000 2.4800000 0.3800000 0.7400000
## Strawberry 0.5714286  1.1904762    0.6666667 1.6666667 0.5714286 0.4285714
##                Grape      Kiwi Raspberry  Rosehips Strawberry    Tomato
## Blackberry 0.2068966 0.0862069 0.4655172 0.3793103  0.5862069 0.2241379
## Cherry     0.9700000 1.0500000 1.5800000 1.2400000  0.6700000 0.6969697
## Strawberry 0.4761905 0.1428571 0.6666667 0.3809524  0.8571429 0.0952381
## 
## , , G2
## 
##              Apricot Blackberry Blackcurrant   Cherry Cranberry      Fig
## Blackberry  6.826531   24.62245     18.61224 23.41837  8.428571 12.73469
## Cherry     10.423077   16.09615     19.09615 30.38462  7.480769 13.46154
## Strawberry 15.450000   21.92500     22.12500 24.92500  7.200000 13.02500
##               Grape     Kiwi Raspberry Rosehips Strawberry   Tomato
## Blackberry 24.55102 16.57143  16.59184 11.23469   6.938776 12.53061
## Cherry     11.67308 11.90385  12.78846 11.46154  10.673077 10.59615
## Strawberry 12.55000 14.87500  19.45000 16.20000  12.775000 14.90000
ggplot2::ggplot(data = data_PREF, 
                aes(x = Test_environment, y = Nb_eggs, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PREF, 
                aes(x = Nb_eggs, fill = Original_environment)) +
  facet_wrap( ~ Generation) +
  geom_histogram(position="identity", alpha=0.5) +
  theme_LO_sober

5.1.2 Analysis

#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment,
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PREF)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[4,2])/(1/anova(m1)[4, 1]))
## [1] 10.1742
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[4, 1])) 
## [1] 0.004407544
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            BoxID, 
          data = data_PREF)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 9.833047
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[6, 1]) )
## [1] 0.004993505
## F test for SA
(Fratio_NonGen <- (anova(m2)[5,2]/anova(m2)[7,2])/(1/anova(m2)[7, 1]))
## [1] 2.500166
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[7, 1]) )
## [1] 0.1287796
## Compute R2 for SA 
## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
(r2_SA_genet <- 1-(anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))))
## [1] 0.2864799
(r2_SA_nongenet <- 1-(anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))))
## [1] 0.0638364
#Local adaptation pattern: 
lm_val = lm(log(Nb_eggs+1) ~ Test_environment + Population + SA + 
                      Test_environment:Original_environment + BoxID, 
                    data = data_PREF[data_PREF$Generation=="G0",])

(Fratio = anova(lm_val)[3,3]/anova(lm_val)[5,3])
## [1] 4.476081
(pvalue = 1 - pf(Fratio,anova(lm_val)[3,1],anova(lm_val)[5,1]))
## [1] 0.0464917
(df1 = anova(lm_val)[3,1])
## [1] 1
(df2 = anova(lm_val)[5,1])
## [1] 21
lm_val = lm(log(Nb_eggs+1) ~ Test_environment + Population + SA + 
                      Test_environment:Original_environment + BoxID, 
                    data = data_PREF[data_PREF$Generation=="G2",])

(Fratio = anova(lm_val)[3,3]/anova(lm_val)[5,3])
## [1] 6.319972
(pvalue = 1 - pf(Fratio,anova(lm_val)[3,1],anova(lm_val)[5,1]))
## [1] 0.02016054
(df1 = anova(lm_val)[3,1])
## [1] 1
(df2 = anova(lm_val)[5,1])
## [1] 21

5.2 Three fruits

5.2.1 Exploration data

tapply(data_PREF_three$Nb_eggs, list(data_PREF_three$Original_environment, 
                                 data_PREF_three$Test_environment, 
                                 data_PREF_three$Generation), mean, na.rm = TRUE)
## , , G0
## 
##            Blackberry    Cherry Strawberry
## Blackberry  0.7931034 0.4137931  0.5862069
## Cherry      0.5000000 2.4800000  0.6700000
## Strawberry  1.1904762 1.6666667  0.8571429
## 
## , , G2
## 
##            Blackberry   Cherry Strawberry
## Blackberry   24.62245 23.41837   6.938776
## Cherry       16.09615 30.38462  10.673077
## Strawberry   21.92500 24.92500  12.775000
ggplot2::ggplot(data = data_PREF_three, 
                aes(x = Test_environment, y = Nb_eggs, color = Test_environment)) +
  facet_wrap( ~ Population) +
  geom_point() +
  geom_boxplot() +
  theme_LO_sober

ggplot2::ggplot(data = data_PREF_three, 
                aes(x = Nb_eggs, fill = Original_environment)) +
  facet_wrap( ~ Generation) +
  geom_histogram(position="identity", alpha=0.5) +
  theme_LO_sober

5.2.2 Analysis

#######################################################
## Analysis of genetic effects lm   ###
#######################################################
  
m1 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA + 
            Original_environment:Test_environment +
            BoxID,
          contrasts = list(Original_environment = "contr.sum", 
                           Test_environment = "contr.sum"), 
          data = data_PREF_three)
  
  
## F test for SA
(Fratio <- (anova(m1)[3,2]/anova(m1)[5,2])/(1/anova(m1)[5, 1]))
## [1] 18.38385
(pvalue <- 1 - pf(Fratio, 1, anova(m1)[5, 1])) 
## [1] 0.02331823
#######################################################
## Analysis of non-genetic effects lm   ###
#######################################################
  
m2 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            BoxID, 
          data = data_PREF_three)

## F test for SA
(Fratio_Gen <- (anova(m2)[3,2]/anova(m2)[6,2])/(1/anova(m2)[6, 1]))
## [1] 12.77418
(pvalue_Gen <- 1 - pf(Fratio_Gen, 1, anova(m2)[6, 1]) )
## [1] 0.0374429
## F test for SA
(Fratio_NonGen <- (anova(m2)[5,2]/anova(m2)[7,2])/(1/anova(m2)[7, 1]))
## [1] 3.536493
(pvalue_NonGen <- 1 - pf(Fratio_NonGen, 1, anova(m2)[7, 1]) )
## [1] 0.1566089
## Compute R2 for SA 
# ## = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
# (r2_SA_genet <- 1-(anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))))
# 
# (r2_SA_nongenet <- 1-(anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))))
  ## Compute R2 = MS Interaction model without SA - MS Interaction model with SA / MS Interaction model without SA
rsqgen <- 1-anova(m2)[6, 3]/((anova(m2)[3, 2]+anova(m2)[6, 2])/(anova(m2)[3, 1]+anova(m2)[6, 1]))
rsqng <- 1-anova(m2)[7, 3]/((anova(m2)[5, 2]+anova(m2)[7, 2])/(anova(m2)[5, 1]+anova(m2)[7, 1]))
  rsqgen
## [1] 0.746421
  rsqng
## [1] 0.3880511
#######################################################
## Extract SA value   ###
#######################################################
#Model used previously
m2 <- aov(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA +
            Original_environment:Test_environment + 
            Original_environment:Test_environment:IndicG0 +
            BoxID, 
          data = data_PREF_three)
summary(m2)
##                                                Df Sum Sq Mean Sq F value
## pop_gen                                        46 1243.4   27.03  36.767
## hab_gen                                         4  135.1   33.78  45.940
## SA                                              1   14.8   14.82  20.162
## BoxID                                         322  405.3    1.26   1.712
## SA:IndicG0                                      1    3.9    3.85   5.237
## Original_environment:Test_environment           3    3.5    1.16   1.578
## IndicG0:Original_environment:Test_environment   3    3.3    1.09   1.481
## Residuals                                     726  533.8    0.74        
##                                                 Pr(>F)    
## pop_gen                                        < 2e-16 ***
## hab_gen                                        < 2e-16 ***
## SA                                            8.28e-06 ***
## BoxID                                         2.34e-09 ***
## SA:IndicG0                                      0.0224 *  
## Original_environment:Test_environment           0.1933    
## IndicG0:Original_environment:Test_environment   0.2184    
## Residuals                                                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#but without interaction to have correct estimate of SA 
m3 <- lm(log(Nb_eggs+1) ~ pop_gen + hab_gen +  SA:IndicG0 + SA +
           BoxID, 
          data = data_PREF_three)
# Pour les SA:
# summary(m3)
# coef(m3)
cf <-coef(summary(m3,complete = TRUE)) 
indic <- grep("SA",rownames(cf))
SAcoef <- cf[indic,1]
names(SAcoef) <- c("SAGen", "SANonGen")
SAcoef
##     SAGen  SANonGen 
## 0.3842848 0.2660122
#With Nico
m3 <- lm(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA:IndicG0 + SA + BoxID,
data = data_PREF_three)
summary(m3)
## 
## Call:
## lm(formula = log(Nb_eggs + 1) ~ pop_gen + hab_gen + SA:IndicG0 + 
##     SA + BoxID, data = data_PREF_three)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.7761 -0.3517 -0.0258  0.3229  3.0050 
## 
## Coefficients: (48 not defined because of singularities)
##                          Estimate Std. Error t value Pr(>|t|)    
## (Intercept)            -3.585e-01  5.084e-01  -0.705 0.481001    
## pop_genBlackberry31_G2  4.903e-01  7.105e-01   0.690 0.490366    
## pop_genBlackberry32_G0  2.310e-01  7.016e-01   0.329 0.742017    
## pop_genBlackberry32_G2  1.954e+00  7.105e-01   2.750 0.006101 ** 
## pop_genBlackberry33_G0  3.662e-01  7.016e-01   0.522 0.601869    
## pop_genBlackberry33_G2  1.888e+00  7.105e-01   2.657 0.008064 ** 
## pop_genBlackberry34_G0  2.855e-14  7.016e-01   0.000 1.000000    
## pop_genBlackberry34_G2  1.494e+00  7.105e-01   2.103 0.035840 *  
## pop_genBlackberry35_G0  7.225e-14  7.016e-01   0.000 1.000000    
## pop_genBlackberry35_G2  2.358e+00  7.105e-01   3.319 0.000949 ***
## pop_genBlackberry36_G0  9.635e-01  7.016e-01   1.373 0.170113    
## pop_genBlackberry36_G2  3.060e+00  7.105e-01   4.307 1.88e-05 ***
## pop_genBlackberry37_G0  2.282e-14  7.016e-01   0.000 1.000000    
## pop_genBlackberry37_G2  2.500e+00  7.105e-01   3.519 0.000461 ***
## pop_genBlackberry38_G0  3.821e-14  7.016e-01   0.000 1.000000    
## pop_genBlackberry38_G2  2.787e+00  7.105e-01   3.922 9.60e-05 ***
## pop_genBlackberry39_G0  6.931e-01  7.016e-01   0.988 0.323515    
## pop_genBlackberry39_G2  5.380e-01  7.105e-01   0.757 0.449169    
## pop_genBlackberry40_G0  1.426e+00  7.016e-01   2.032 0.042534 *  
## pop_genBlackberry40_G2  1.075e+00  7.105e-01   1.512 0.130904    
## pop_genBlackberry43_G2  1.992e+00  7.105e-01   2.803 0.005189 ** 
## pop_genBlackberry44_G0  2.310e-01  7.016e-01   0.329 0.742017    
## pop_genBlackberry44_G2  3.094e+00  7.105e-01   4.355 1.52e-05 ***
## pop_genBlackberry45_G0  4.858e-14  7.016e-01   0.000 1.000000    
## pop_genBlackberry45_G2  1.435e+00  7.105e-01   2.019 0.043829 *  
## pop_genCherry103_G0     3.662e-01  7.016e-01   0.522 0.601869    
## pop_genCherry103_G2     9.771e-01  7.105e-01   1.375 0.169493    
## pop_genCherry104_G0     5.365e-01  7.016e-01   0.765 0.444738    
## pop_genCherry104_G2     6.502e-01  7.105e-01   0.915 0.360463    
## pop_genCherry3_G0       1.397e+00  7.016e-01   1.990 0.046911 *  
## pop_genCherry3_G2       2.620e+00  7.105e-01   3.687 0.000244 ***
## pop_genCherry47_G0      3.609e-14  7.016e-01   0.000 1.000000    
## pop_genCherry47_G2      9.744e-02  7.105e-01   0.137 0.890963    
## pop_genCherry50_G0      2.310e-01  7.016e-01   0.329 0.742017    
## pop_genCherry50_G2      3.247e+00  7.105e-01   4.570 5.72e-06 ***
## pop_genCherry51_G0      3.166e-14  7.016e-01   0.000 1.000000    
## pop_genCherry52_G0      2.812e-14  7.016e-01   0.000 1.000000    
## pop_genCherry52_G2      1.558e+00  7.105e-01   2.193 0.028631 *  
## pop_genCherry6_G0       3.672e-14  7.016e-01   0.000 1.000000    
## pop_genCherry6_G2       4.030e+00  7.105e-01   5.671 2.04e-08 ***
## pop_genCherry7_G2       2.138e+00  7.105e-01   3.009 0.002708 ** 
## pop_genStrawberry42_G0  7.324e-01  7.016e-01   1.044 0.296884    
## pop_genStrawberry42_G2  1.096e+00  7.105e-01   1.543 0.123380    
## pop_genStrawberry44_G0  2.310e-01  7.016e-01   0.329 0.742017    
## pop_genStrawberry44_G2  2.623e+00  7.105e-01   3.692 0.000239 ***
## pop_genStrawberry53_G0  1.134e+00  7.016e-01   1.616 0.106551    
## pop_genStrawberry53_G2  1.653e+00  7.105e-01   2.327 0.020232 *  
## hab_genBlackberry_G2    6.558e-01  9.144e-02   7.172 1.81e-12 ***
## hab_genCherry_G0        1.452e-01  9.300e-02   1.561 0.118901    
## hab_genCherry_G2        1.129e+00  8.831e-02  12.788  < 2e-16 ***
## hab_genStrawberry_G0    1.394e-02  9.252e-02   0.151 0.880275    
## hab_genStrawberry_G2           NA         NA      NA       NA    
## SA1                     3.843e-01  7.951e-02   4.833 1.64e-06 ***
## BoxID557               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID558                5.825e-15  7.016e-01   0.000 1.000000    
## BoxID559               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID560               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID561               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID562                1.352e-01  7.016e-01   0.193 0.847300    
## BoxID563               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID564                4.950e-15  7.016e-01   0.000 1.000000    
## BoxID565                1.086e+00  7.016e-01   1.548 0.122080    
## BoxID566               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID567                6.716e-01  7.016e-01   0.957 0.338750    
## BoxID568               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID577                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID578                1.484e-15  7.016e-01   0.000 1.000000    
## BoxID579                1.122e-15  7.016e-01   0.000 1.000000    
## BoxID580                1.076e-15  7.016e-01   0.000 1.000000    
## BoxID581                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID582                1.127e-15  7.016e-01   0.000 1.000000    
## BoxID583                4.426e-16  7.016e-01   0.000 1.000000    
## BoxID584                1.627e-15  7.016e-01   0.000 1.000000    
## BoxID585                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID586                1.715e-15  7.016e-01   0.000 1.000000    
## BoxID587                1.229e-15  7.016e-01   0.000 1.000000    
## BoxID588                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID589                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID590                1.646e-15  7.016e-01   0.000 1.000000    
## BoxID591                1.538e-15  7.016e-01   0.000 1.000000    
## BoxID592                       NA         NA      NA       NA    
## BoxID593                5.982e-15  7.016e-01   0.000 1.000000    
## BoxID594                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID595                5.890e-15  7.016e-01   0.000 1.000000    
## BoxID596                5.887e-15  7.016e-01   0.000 1.000000    
## BoxID597                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID598                5.841e-15  7.016e-01   0.000 1.000000    
## BoxID599                6.204e-15  7.016e-01   0.000 1.000000    
## BoxID600                5.973e-01  7.016e-01   0.851 0.394909    
## BoxID604                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID605                6.241e-15  7.016e-01   0.000 1.000000    
## BoxID606                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID607                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID608                6.024e-15  7.016e-01   0.000 1.000000    
## BoxID609                7.325e-15  7.016e-01   0.000 1.000000    
## BoxID610                       NA         NA      NA       NA    
## BoxID611               -1.423e-15  7.016e-01   0.000 1.000000    
## BoxID612                3.662e-01  7.016e-01   0.522 0.601869    
## BoxID613                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID614                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID615               -1.012e-15  7.016e-01   0.000 1.000000    
## BoxID616               -1.378e-15  7.016e-01   0.000 1.000000    
## BoxID617                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID620               -4.939e-01  7.016e-01  -0.704 0.481720    
## BoxID621               -9.345e-01  7.016e-01  -1.332 0.183323    
## BoxID622               -1.166e+00  7.016e-01  -1.661 0.097109 .  
## BoxID623               -1.166e+00  7.016e-01  -1.661 0.097109 .  
## BoxID624               -1.397e+00  7.016e-01  -1.990 0.046911 *  
## BoxID625                       NA         NA      NA       NA    
## BoxID629                4.849e-15  7.016e-01   0.000 1.000000    
## BoxID630                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID631                4.377e-15  7.016e-01   0.000 1.000000    
## BoxID632                3.940e-15  7.016e-01   0.000 1.000000    
## BoxID633                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID634                3.666e-15  7.016e-01   0.000 1.000000    
## BoxID635                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID636                3.662e-01  7.016e-01   0.522 0.601869    
## BoxID638                4.057e-15  7.016e-01   0.000 1.000000    
## BoxID639                3.439e-15  7.016e-01   0.000 1.000000    
## BoxID640                3.435e-15  7.016e-01   0.000 1.000000    
## BoxID641                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID642                       NA         NA      NA       NA    
## BoxID643                8.446e-01  7.016e-01   1.204 0.229080    
## BoxID644                3.877e-01  7.016e-01   0.553 0.580704    
## BoxID650                       NA         NA      NA       NA    
## BoxID653               -1.124e-15  7.016e-01   0.000 1.000000    
## BoxID654               -1.242e-15  7.016e-01   0.000 1.000000    
## BoxID655                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID656               -1.198e-15  7.016e-01   0.000 1.000000    
## BoxID657               -1.196e-15  7.016e-01   0.000 1.000000    
## BoxID658               -1.226e-15  7.016e-01   0.000 1.000000    
## BoxID660               -1.374e-15  7.016e-01   0.000 1.000000    
## BoxID661                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID662                       NA         NA      NA       NA    
## BoxID665                4.184e-16  7.016e-01   0.000 1.000000    
## BoxID666               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID667               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID668               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID669               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID670               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID671                8.338e-16  7.016e-01   0.000 1.000000    
## BoxID673               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID674                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID675               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID676               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID677                       NA         NA      NA       NA    
## BoxID678                1.147e+00  7.016e-01   1.635 0.102581    
## BoxID679                6.077e-02  7.016e-01   0.087 0.930998    
## BoxID680                1.297e+00  7.016e-01   1.849 0.064865 .  
## BoxID682                6.486e-01  7.016e-01   0.924 0.355538    
## BoxID683                6.077e-02  7.016e-01   0.087 0.930998    
## BoxID684                6.077e-02  7.016e-01   0.087 0.930998    
## BoxID685                9.850e-01  7.016e-01   1.404 0.160786    
## BoxID686                5.743e-01  7.016e-01   0.818 0.413354    
## BoxID687                4.621e-01  7.016e-01   0.659 0.510349    
## BoxID688               -3.054e-01  7.016e-01  -0.435 0.663457    
## BoxID689                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID690               -3.054e-01  7.016e-01  -0.435 0.663457    
## BoxID691                1.547e+00  7.016e-01   2.205 0.027787 *  
## BoxID692                3.185e-01  7.016e-01   0.454 0.649995    
## BoxID693                3.662e-01  7.016e-01   0.522 0.601869    
## BoxID694                5.087e-01  7.016e-01   0.725 0.468673    
## BoxID695                1.018e+00  7.016e-01   1.451 0.147227    
## BoxID696                       NA         NA      NA       NA    
## BoxID775               -7.324e-01  7.016e-01  -1.044 0.296884    
## BoxID776               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID777                1.703e-01  7.016e-01   0.243 0.808314    
## BoxID778                       NA         NA      NA       NA    
## BoxID779               -7.324e-01  7.016e-01  -1.044 0.296884    
## BoxID780               -1.352e-01  7.016e-01  -0.193 0.847300    
## BoxID782               -5.014e-01  7.016e-01  -0.715 0.475100    
## BoxID783                3.269e-01  7.016e-01   0.466 0.641366    
## BoxID784                       NA         NA      NA       NA    
## BoxID794                3.662e-01  7.016e-01   0.522 0.601869    
## BoxID795               -1.348e+00  7.016e-01  -1.921 0.055171 .  
## BoxID796               -1.776e+00  7.016e-01  -2.531 0.011590 *  
## BoxID797               -3.937e-01  7.016e-01  -0.561 0.574842    
## BoxID798                1.477e+00  7.016e-01   2.105 0.035628 *  
## BoxID800               -9.413e-01  7.016e-01  -1.342 0.180141    
## BoxID806                       NA         NA      NA       NA    
## BoxID814               -7.473e-01  7.016e-01  -1.065 0.287180    
## BoxID815               -2.965e-01  7.016e-01  -0.423 0.672716    
## BoxID816               -8.298e-01  7.016e-01  -1.183 0.237314    
## BoxID817                5.391e-03  7.016e-01   0.008 0.993872    
## BoxID818               -4.362e-01  7.016e-01  -0.622 0.534316    
## BoxID819                       NA         NA      NA       NA    
## BoxID820                       NA         NA      NA       NA    
## BoxID828               -1.623e+00  7.016e-01  -2.314 0.020968 *  
## BoxID829               -1.003e-01  7.016e-01  -0.143 0.886355    
## BoxID831               -2.503e+00  7.016e-01  -3.567 0.000384 ***
## BoxID832                9.842e-01  7.016e-01   1.403 0.161097    
## BoxID833               -1.430e+00  7.016e-01  -2.038 0.041899 *  
## BoxID834               -1.623e+00  7.016e-01  -2.314 0.020968 *  
## BoxID835               -8.463e-01  7.016e-01  -1.206 0.228112    
## BoxID836                2.735e-02  7.016e-01   0.039 0.968910    
## BoxID837                6.586e-01  7.016e-01   0.939 0.348190    
## BoxID838               -1.273e+00  7.016e-01  -1.815 0.069963 .  
## BoxID839               -4.548e-01  7.016e-01  -0.648 0.517007    
## BoxID840               -1.202e-01  7.016e-01  -0.171 0.863983    
## BoxID841               -5.599e-01  7.016e-01  -0.798 0.425145    
## BoxID842               -9.302e-01  7.016e-01  -1.326 0.185314    
## BoxID843               -1.415e+00  7.016e-01  -2.017 0.044038 *  
## BoxID851                6.393e-01  7.016e-01   0.911 0.362470    
## BoxID852               -1.003e+00  7.016e-01  -1.430 0.153147    
## BoxID853                4.544e-01  7.016e-01   0.648 0.517436    
## BoxID854                1.180e+00  7.016e-01   1.682 0.092922 .  
## BoxID855               -8.629e-01  7.016e-01  -1.230 0.219140    
## BoxID856               -3.855e-01  7.016e-01  -0.549 0.582836    
## BoxID857                6.130e-02  7.016e-01   0.087 0.930407    
## BoxID858               -9.588e-01  7.016e-01  -1.367 0.172187    
## BoxID859               -4.765e-01  7.016e-01  -0.679 0.497272    
## BoxID860               -1.352e-01  7.016e-01  -0.193 0.847300    
## BoxID861                1.665e+00  7.016e-01   2.372 0.017926 *  
## BoxID862                1.003e+00  7.016e-01   1.429 0.153388    
## BoxID863                1.451e+00  7.016e-01   2.068 0.039033 *  
## BoxID864                2.060e+00  7.016e-01   2.936 0.003431 ** 
## BoxID865                2.271e+00  7.016e-01   3.237 0.001261 ** 
## BoxID866                2.616e+00  7.016e-01   3.728 0.000208 ***
## BoxID867                2.007e+00  7.016e-01   2.861 0.004342 ** 
## BoxID869                       NA         NA      NA       NA    
## BoxID870                4.039e-01  7.016e-01   0.576 0.565019    
## BoxID871               -9.588e-01  7.016e-01  -1.367 0.172187    
## BoxID872                6.646e-02  7.016e-01   0.095 0.924557    
## BoxID873                3.518e-01  7.016e-01   0.501 0.616219    
## BoxID874                5.286e-01  7.016e-01   0.753 0.451485    
## BoxID875               -5.412e-01  7.016e-01  -0.771 0.440737    
## BoxID876               -7.493e-01  7.016e-01  -1.068 0.285917    
## BoxID877               -1.421e+00  7.016e-01  -2.025 0.043213 *  
## BoxID878               -4.728e-03  7.016e-01  -0.007 0.994625    
## BoxID881               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID882               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID883                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID884               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID885               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID886               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID887               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID888               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID889               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID890                7.203e-01  7.016e-01   1.027 0.304945    
## BoxID891               -6.486e-01  7.016e-01  -0.924 0.355538    
## BoxID893               -8.797e-01  7.016e-01  -1.254 0.210317    
## BoxID894                1.565e+00  7.016e-01   2.230 0.026046 *  
## BoxID895                8.253e-01  7.016e-01   1.176 0.239858    
## BoxID902               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID904               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID905                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID906               -2.310e-01  7.016e-01  -0.329 0.742017    
## BoxID907                2.609e-14  7.016e-01   0.000 1.000000    
## BoxID908                       NA         NA      NA       NA    
## BoxID910                2.910e-01  7.016e-01   0.415 0.678475    
## BoxID911                3.420e+00  7.016e-01   4.874 1.34e-06 ***
## BoxID912                2.694e+00  7.016e-01   3.840 0.000134 ***
## BoxID913               -5.966e-01  7.016e-01  -0.850 0.395409    
## BoxID914                       NA         NA      NA       NA    
## BoxID915               -1.043e+00  7.016e-01  -1.486 0.137611    
## BoxID916                2.280e+00  7.016e-01   3.250 0.001207 ** 
## BoxID917                1.301e+00  7.016e-01   1.854 0.064171 .  
## BoxID918                1.749e+00  7.016e-01   2.493 0.012893 *  
## BoxID919                2.640e+00  7.016e-01   3.763 0.000181 ***
## BoxID920                2.751e+00  7.016e-01   3.921 9.64e-05 ***
## BoxID921                1.551e+00  7.016e-01   2.211 0.027340 *  
## BoxID922                1.166e+00  7.016e-01   1.661 0.097109 .  
## BoxID923                6.826e-01  7.016e-01   0.973 0.330953    
## BoxID924                9.999e-01  7.016e-01   1.425 0.154526    
## BoxID925                       NA         NA      NA       NA    
## BoxID926                3.128e-01  7.016e-01   0.446 0.655850    
## BoxID927               -1.735e+00  7.016e-01  -2.473 0.013608 *  
## BoxID928                       NA         NA      NA       NA    
## BoxID933                1.590e+00  7.016e-01   2.266 0.023732 *  
## BoxID934                7.336e-01  7.016e-01   1.046 0.296084    
## BoxID935                1.043e+00  7.016e-01   1.486 0.137646    
## BoxID936                1.260e+00  7.016e-01   1.796 0.072893 .  
## BoxID937                       NA         NA      NA       NA    
## BoxID939                       NA         NA      NA       NA    
## BoxID943                1.516e+00  7.016e-01   2.161 0.031028 *  
## BoxID944               -1.865e-01  7.016e-01  -0.266 0.790415    
## BoxID945               -3.662e-01  7.016e-01  -0.522 0.601869    
## BoxID946                1.022e+00  7.016e-01   1.456 0.145814    
## BoxID947                       NA         NA      NA       NA    
## BoxID949                5.661e-01  7.016e-01   0.807 0.420034    
## BoxID950                       NA         NA      NA       NA    
## BoxID951               -1.326e-14  7.016e-01   0.000 1.000000    
## BoxID952                8.283e-01  7.016e-01   1.181 0.238161    
## BoxID953               -7.637e-14  7.016e-01   0.000 1.000000    
## BoxID954                6.931e-01  7.016e-01   0.988 0.323515    
## BoxID955                6.931e-01  7.016e-01   0.988 0.323515    
## BoxID956                7.416e-15  7.016e-01   0.000 1.000000    
## BoxID957               -5.306e-15  7.016e-01   0.000 1.000000    
## BoxID958                2.172e-14  7.016e-01   0.000 1.000000    
## BoxID959                       NA         NA      NA       NA    
## BoxID960                       NA         NA      NA       NA    
## BoxID962                3.128e-01  7.016e-01   0.446 0.655899    
## BoxID963                       NA         NA      NA       NA    
## BoxID965                1.398e-15  7.016e-01   0.000 1.000000    
## BoxID966                2.581e-14  7.016e-01   0.000 1.000000    
## BoxID967                9.259e-16  7.016e-01   0.000 1.000000    
## BoxID968                1.505e-15  7.016e-01   0.000 1.000000    
## BoxID969                9.872e-16  7.016e-01   0.000 1.000000    
## BoxID970                2.202e-15  7.016e-01   0.000 1.000000    
## BoxID971                1.242e-15  7.016e-01   0.000 1.000000    
## BoxID972                       NA         NA      NA       NA    
## BoxID973                2.680e-14  7.016e-01   0.000 1.000000    
## BoxID974                2.597e-14  7.016e-01   0.000 1.000000    
## BoxID975                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID976                1.230e+00  7.016e-01   1.753 0.080097 .  
## BoxID977                2.310e-01  7.016e-01   0.329 0.742017    
## BoxID978                       NA         NA      NA       NA    
## BoxID979                3.662e-01  7.016e-01   0.522 0.601869    
## BoxID980                7.993e-01  7.016e-01   1.139 0.254984    
## BoxID983               -6.931e-01  7.016e-01  -0.988 0.323515    
## BoxID984               -1.671e-15  7.016e-01   0.000 1.000000    
## BoxID985                       NA         NA      NA       NA    
## BoxID986                5.365e-01  7.016e-01   0.765 0.444738    
## BoxID987                5.634e-15  7.016e-01   0.000 1.000000    
## BoxID988                2.542e-15  7.016e-01   0.000 1.000000    
## BoxID989                       NA         NA      NA       NA    
## BoxID990                1.352e-01  7.016e-01   0.193 0.847300    
## BoxID991                2.613e-14  7.016e-01   0.000 1.000000    
## BoxID992                       NA         NA      NA       NA    
## BoxID993                       NA         NA      NA       NA    
## BoxID994               -1.426e+00  7.016e-01  -2.032 0.042534 *  
## BoxID997                       NA         NA      NA       NA    
## BoxID998                       NA         NA      NA       NA    
## BoxID999                       NA         NA      NA       NA    
## BoxID1001               2.293e+00  7.016e-01   3.268 0.001133 ** 
## BoxID1002               3.105e+00  7.016e-01   4.426 1.11e-05 ***
## BoxID1003               3.105e-01  7.016e-01   0.443 0.658203    
## BoxID1004               4.621e-01  7.016e-01   0.659 0.510349    
## BoxID1005               1.676e+00  7.016e-01   2.389 0.017159 *  
## BoxID1006               2.021e+00  7.016e-01   2.880 0.004093 ** 
## BoxID1008               1.338e+00  7.016e-01   1.907 0.056891 .  
## BoxID1009               1.890e+00  7.016e-01   2.694 0.007220 ** 
## BoxID1010                      NA         NA      NA       NA    
## BoxID1011               1.843e+00  7.016e-01   2.627 0.008795 ** 
## BoxID1012               2.305e+00  7.016e-01   3.286 0.001066 ** 
## BoxID1013               4.784e-01  7.016e-01   0.682 0.495583    
## BoxID1014               1.873e+00  7.016e-01   2.669 0.007775 ** 
## BoxID1015              -8.582e-01  7.016e-01  -1.223 0.221673    
## BoxID1016              -7.068e-01  7.016e-01  -1.007 0.314114    
## BoxID1017              -1.259e+00  7.016e-01  -1.795 0.073045 .  
## BoxID1018              -1.586e+00  7.016e-01  -2.261 0.024045 *  
## BoxID1019              -1.073e+00  7.016e-01  -1.529 0.126632    
## BoxID1020              -2.319e+00  7.016e-01  -3.305 0.000996 ***
## BoxID1021               5.409e-01  7.016e-01   0.771 0.441001    
## BoxID1022                      NA         NA      NA       NA    
## BoxID1024              -1.258e-01  7.016e-01  -0.179 0.857791    
## BoxID1025               4.791e-01  7.016e-01   0.683 0.494903    
## BoxID1026              -1.857e+00  7.016e-01  -2.646 0.008311 ** 
## BoxID1027              -1.626e+00  7.016e-01  -2.317 0.020775 *  
## BoxID1028              -9.431e-01  7.016e-01  -1.344 0.179290    
## BoxID1029              -6.169e-01  7.016e-01  -0.879 0.379579    
## BoxID1030                      NA         NA      NA       NA    
## BoxID1036               1.600e+00  7.016e-01   2.280 0.022899 *  
## BoxID1037              -2.903e-01  7.016e-01  -0.414 0.679198    
## BoxID1038               9.744e-02  7.016e-01   0.139 0.889583    
## BoxID1039               9.219e-02  7.016e-01   0.131 0.895496    
## BoxID1040               5.501e-01  7.016e-01   0.784 0.433227    
## BoxID1041               8.073e-01  7.016e-01   1.151 0.250259    
## BoxID1042              -7.916e-01  7.016e-01  -1.128 0.259564    
## BoxID1043               1.322e+00  7.016e-01   1.884 0.060017 .  
## BoxID1047               1.093e+00  7.016e-01   1.558 0.119629    
## BoxID1048              -1.289e+00  7.016e-01  -1.837 0.066619 .  
## BoxID1049               4.867e-01  7.016e-01   0.694 0.488142    
## BoxID1050                      NA         NA      NA       NA    
## BoxID1051              -8.658e-01  7.016e-01  -1.234 0.217586    
## BoxID1052              -1.544e+00  7.016e-01  -2.201 0.028026 *  
## BoxID1056               3.297e-01  7.016e-01   0.470 0.638544    
## BoxID1057               1.879e-02  7.016e-01   0.027 0.978644    
## BoxID1058                      NA         NA      NA       NA    
## BoxID1059               2.260e-01  7.016e-01   0.322 0.747503    
## BoxID1060              -7.993e-01  7.016e-01  -1.139 0.254984    
## BoxID1061              -2.403e+00  7.016e-01  -3.424 0.000650 ***
## BoxID1062               3.889e-01  7.016e-01   0.554 0.579577    
## BoxID1063              -3.264e-01  7.016e-01  -0.465 0.641869    
## BoxID1064              -1.589e+00  7.016e-01  -2.265 0.023810 *  
## BoxID1065               8.989e-02  7.016e-01   0.128 0.898093    
## BoxID1066              -4.156e-01  7.016e-01  -0.592 0.553794    
## BoxID1067              -1.337e+00  7.016e-01  -1.905 0.057181 .  
## BoxID1068              -3.262e-01  7.016e-01  -0.465 0.642136    
## BoxID1069               3.136e-02  7.016e-01   0.045 0.964364    
## BoxID1070               2.850e-01  7.016e-01   0.406 0.684739    
## BoxID1071                      NA         NA      NA       NA    
## BoxID1072              -9.217e-01  7.016e-01  -1.314 0.189349    
## BoxID1073              -5.917e-01  7.016e-01  -0.843 0.399282    
## BoxID1074              -1.286e-01  7.016e-01  -0.183 0.854672    
## BoxID1075              -5.813e-01  7.016e-01  -0.829 0.407616    
## BoxID1076              -5.014e-01  7.016e-01  -0.715 0.475100    
## BoxID1078              -4.577e-01  7.016e-01  -0.652 0.514397    
## BoxID1079              -2.097e+00  7.016e-01  -2.989 0.002892 ** 
## BoxID1080                      NA         NA      NA       NA    
## BoxID1081              -2.840e-01  7.016e-01  -0.405 0.685780    
## BoxID1085               6.433e-01  7.016e-01   0.917 0.359505    
## BoxID1086               2.721e-01  7.016e-01   0.388 0.698296    
## BoxID1087              -3.946e-01  7.016e-01  -0.562 0.574016    
## BoxID1088              -8.458e-01  7.016e-01  -1.205 0.228424    
## BoxID1089               8.222e-02  7.016e-01   0.117 0.906739    
## BoxID1090                      NA         NA      NA       NA    
## BoxID1091              -2.467e-01  7.016e-01  -0.352 0.725209    
## BoxID1092              -5.527e-01  7.016e-01  -0.788 0.431063    
## BoxID1093                      NA         NA      NA       NA    
## BoxID1094                      NA         NA      NA       NA    
## BoxID1095              -8.258e-01  7.016e-01  -1.177 0.239594    
## BoxID1096              -7.406e-01  7.016e-01  -1.056 0.291540    
## BoxID1097              -1.886e+00  7.016e-01  -2.688 0.007350 ** 
## BoxID1098              -7.403e-02  7.016e-01  -0.106 0.916000    
## BoxID1101              -2.399e+00  7.016e-01  -3.420 0.000661 ***
## BoxID1102              -1.776e+00  7.016e-01  -2.531 0.011594 *  
## BoxID1103              -2.992e-01  7.016e-01  -0.427 0.669868    
## BoxID1104              -4.735e-02  7.016e-01  -0.067 0.946216    
## BoxID1105              -1.003e+00  7.016e-01  -1.429 0.153303    
## BoxID1107              -3.659e-01  7.016e-01  -0.521 0.602219    
## BoxID1108               5.870e-01  7.016e-01   0.837 0.403045    
## BoxID1109              -1.371e+00  7.016e-01  -1.954 0.051070 .  
## BoxID1110                      NA         NA      NA       NA    
## BoxID1112               2.114e+00  7.016e-01   3.013 0.002672 ** 
## BoxID1115              -2.300e-02  7.016e-01  -0.033 0.973861    
## BoxID1116               1.200e+00  7.016e-01   1.710 0.087656 .  
## BoxID1117              -1.778e-01  7.016e-01  -0.253 0.800057    
## BoxID1118              -7.438e-02  7.016e-01  -0.106 0.915601    
## BoxID1119               1.729e-01  7.016e-01   0.246 0.805384    
## BoxID1120                      NA         NA      NA       NA    
## SA0:IndicG0             2.660e-01  1.165e-01   2.284 0.022683 *  
## SA1:IndicG0                    NA         NA      NA       NA    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8593 on 732 degrees of freedom
## Multiple R-squared:  0.7693, Adjusted R-squared:  0.6514 
## F-statistic: 6.527 on 374 and 732 DF,  p-value: < 2.2e-16
coef(m3)
##            (Intercept) pop_genBlackberry31_G2 pop_genBlackberry32_G0 
##          -3.584823e-01           4.903231e-01           2.310491e-01 
## pop_genBlackberry32_G2 pop_genBlackberry33_G0 pop_genBlackberry33_G2 
##           1.954188e+00           3.662041e-01           1.887651e+00 
## pop_genBlackberry34_G0 pop_genBlackberry34_G2 pop_genBlackberry35_G0 
##           2.855444e-14           1.493990e+00           7.224669e-14 
## pop_genBlackberry35_G2 pop_genBlackberry36_G0 pop_genBlackberry36_G2 
##           2.358169e+00           9.634573e-01           3.060449e+00 
## pop_genBlackberry37_G0 pop_genBlackberry37_G2 pop_genBlackberry38_G0 
##           2.281848e-14           2.500058e+00           3.821188e-14 
## pop_genBlackberry38_G2 pop_genBlackberry39_G0 pop_genBlackberry39_G2 
##           2.786906e+00           6.931472e-01           5.380234e-01 
## pop_genBlackberry40_G0 pop_genBlackberry40_G2 pop_genBlackberry43_G2 
##           1.425555e+00           1.074503e+00           1.991965e+00 
## pop_genBlackberry44_G0 pop_genBlackberry44_G2 pop_genBlackberry45_G0 
##           2.310491e-01           3.094169e+00           4.857505e-14 
## pop_genBlackberry45_G2    pop_genCherry103_G0    pop_genCherry103_G2 
##           1.434728e+00           3.662041e-01           9.771239e-01 
##    pop_genCherry104_G0    pop_genCherry104_G2      pop_genCherry3_G0 
##           5.364793e-01           6.501808e-01           1.396552e+00 
##      pop_genCherry3_G2     pop_genCherry47_G0     pop_genCherry47_G2 
##           2.619542e+00           3.609108e-14           9.743814e-02 
##     pop_genCherry50_G0     pop_genCherry50_G2     pop_genCherry51_G0 
##           2.310491e-01           3.247355e+00           3.165834e-14 
##     pop_genCherry52_G0     pop_genCherry52_G2      pop_genCherry6_G0 
##           2.811518e-14           1.558114e+00           3.672278e-14 
##      pop_genCherry6_G2      pop_genCherry7_G2 pop_genStrawberry42_G0 
##           4.029648e+00           2.138269e+00           7.324082e-01 
## pop_genStrawberry42_G2 pop_genStrawberry44_G0 pop_genStrawberry44_G2 
##           1.096016e+00           2.310491e-01           2.623487e+00 
## pop_genStrawberry53_G0 pop_genStrawberry53_G2   hab_genBlackberry_G2 
##           1.133732e+00           1.653486e+00           6.558473e-01 
##       hab_genCherry_G0       hab_genCherry_G2   hab_genStrawberry_G0 
##           1.451974e-01           1.129295e+00           1.394028e-02 
##   hab_genStrawberry_G2                    SA1               BoxID557 
##                     NA           3.842848e-01          -2.310491e-01 
##               BoxID558               BoxID559               BoxID560 
##           5.824560e-15          -2.310491e-01          -2.310491e-01 
##               BoxID561               BoxID562               BoxID563 
##          -2.310491e-01           1.351550e-01          -2.310491e-01 
##               BoxID564               BoxID565               BoxID566 
##           4.950275e-15           1.086032e+00          -2.310491e-01 
##               BoxID567               BoxID568               BoxID577 
##           6.716343e-01          -2.310491e-01           2.310491e-01 
##               BoxID578               BoxID579               BoxID580 
##           1.483963e-15           1.122104e-15           1.075797e-15 
##               BoxID581               BoxID582               BoxID583 
##           4.620981e-01           1.127183e-15           4.426449e-16 
##               BoxID584               BoxID585               BoxID586 
##           1.627439e-15           2.310491e-01           1.714978e-15 
##               BoxID587               BoxID588               BoxID589 
##           1.228792e-15           2.310491e-01           4.620981e-01 
##               BoxID590               BoxID591               BoxID592 
##           1.645725e-15           1.538026e-15                     NA 
##               BoxID593               BoxID594               BoxID595 
##           5.981733e-15           2.310491e-01           5.890445e-15 
##               BoxID596               BoxID597               BoxID598 
##           5.887325e-15           2.310491e-01           5.841222e-15 
##               BoxID599               BoxID600               BoxID604 
##           6.203814e-15           5.972532e-01           2.310491e-01 
##               BoxID605               BoxID606               BoxID607 
##           6.240506e-15           2.310491e-01           2.310491e-01 
##               BoxID608               BoxID609               BoxID610 
##           6.023590e-15           7.325344e-15                     NA 
##               BoxID611               BoxID612               BoxID613 
##          -1.423320e-15           3.662041e-01           2.310491e-01 
##               BoxID614               BoxID615               BoxID616 
##           4.620981e-01          -1.012361e-15          -1.378161e-15 
##               BoxID617               BoxID620               BoxID621 
##           4.620981e-01          -4.938682e-01          -9.344535e-01 
##               BoxID622               BoxID623               BoxID624 
##          -1.165503e+00          -1.165503e+00          -1.396552e+00 
##               BoxID625               BoxID629               BoxID630 
##                     NA           4.849427e-15           2.310491e-01 
##               BoxID631               BoxID632               BoxID633 
##           4.377214e-15           3.939944e-15           2.310491e-01 
##               BoxID634               BoxID635               BoxID636 
##           3.666058e-15           2.310491e-01           3.662041e-01 
##               BoxID638               BoxID639               BoxID640 
##           4.056532e-15           3.439320e-15           3.434577e-15 
##               BoxID641               BoxID642               BoxID643 
##           2.310491e-01                     NA           8.445656e-01 
##               BoxID644               BoxID650               BoxID653 
##           3.877169e-01                     NA          -1.124065e-15 
##               BoxID654               BoxID655               BoxID656 
##          -1.242006e-15           4.620981e-01          -1.198062e-15 
##               BoxID657               BoxID658               BoxID660 
##          -1.196362e-15          -1.225848e-15          -1.373831e-15 
##               BoxID661               BoxID662               BoxID665 
##           2.310491e-01                     NA           4.184310e-16 
##               BoxID666               BoxID667               BoxID668 
##          -2.310491e-01          -2.310491e-01          -2.310491e-01 
##               BoxID669               BoxID670               BoxID671 
##          -2.310491e-01          -2.310491e-01           8.338229e-16 
##               BoxID673               BoxID674               BoxID675 
##          -2.310491e-01           4.620981e-01          -2.310491e-01 
##               BoxID676               BoxID677               BoxID678 
##          -2.310491e-01                     NA           1.146806e+00 
##               BoxID679               BoxID680               BoxID682 
##           6.077385e-02           1.297273e+00           6.486367e-01 
##               BoxID683               BoxID684               BoxID685 
##           6.077385e-02           6.077385e-02           9.849701e-01 
##               BoxID686               BoxID687               BoxID688 
##           5.742555e-01           4.620981e-01          -3.054302e-01 
##               BoxID689               BoxID690               BoxID691 
##           2.310491e-01          -3.054302e-01           1.546846e+00 
##               BoxID692               BoxID693               BoxID694 
##           3.185038e-01           3.662041e-01           5.086854e-01 
##               BoxID695               BoxID696               BoxID775 
##           1.018000e+00                     NA          -7.324082e-01 
##               BoxID776               BoxID777               BoxID778 
##          -3.662041e-01           1.702752e-01                     NA 
##               BoxID779               BoxID780               BoxID782 
##          -7.324082e-01          -1.351550e-01          -5.013591e-01 
##               BoxID783               BoxID784               BoxID794 
##           3.269431e-01                     NA           3.662041e-01 
##               BoxID795               BoxID796               BoxID797 
##          -1.347526e+00          -1.775618e+00          -3.937386e-01 
##               BoxID798               BoxID800               BoxID806 
##           1.476939e+00          -9.412940e-01                     NA 
##               BoxID814               BoxID815               BoxID816 
##          -7.472970e-01          -2.964994e-01          -8.298011e-01 
##               BoxID817               BoxID818               BoxID819 
##           5.390953e-03          -4.362145e-01                     NA 
##               BoxID820               BoxID828               BoxID829 
##                     NA          -1.623243e+00          -1.003091e-01 
##               BoxID831               BoxID832               BoxID833 
##          -2.502929e+00           9.842383e-01          -1.429970e+00 
##               BoxID834               BoxID835               BoxID836 
##          -1.623243e+00          -8.463246e-01           2.735498e-02 
##               BoxID837               BoxID838               BoxID839 
##           6.586153e-01          -1.273303e+00          -4.548460e-01 
##               BoxID840               BoxID841               BoxID842 
##          -1.202340e-01          -5.598730e-01          -9.302179e-01 
##               BoxID843               BoxID851               BoxID852 
##          -1.415314e+00           6.393422e-01          -1.003306e+00 
##               BoxID853               BoxID854               BoxID855 
##           4.543804e-01           1.180388e+00          -8.629011e-01 
##               BoxID856               BoxID857               BoxID858 
##          -3.855332e-01           6.129509e-02          -9.587952e-01 
##               BoxID859               BoxID860               BoxID861 
##          -4.764888e-01          -1.351550e-01           1.664581e+00 
##               BoxID862               BoxID863               BoxID864 
##           1.002718e+00           1.450630e+00           2.059807e+00 
##               BoxID865               BoxID866               BoxID867 
##           2.271377e+00           2.615921e+00           2.007412e+00 
##               BoxID869               BoxID870               BoxID871 
##                     NA           4.038965e-01          -9.587952e-01 
##               BoxID872               BoxID873               BoxID874 
##           6.646316e-02           3.518134e-01           5.285613e-01 
##               BoxID875               BoxID876               BoxID877 
##          -5.412075e-01          -7.492590e-01          -1.420893e+00 
##               BoxID878               BoxID881               BoxID882 
##          -4.728212e-03          -2.310491e-01          -2.310491e-01 
##               BoxID883               BoxID884               BoxID885 
##           2.310491e-01          -3.662041e-01          -3.662041e-01 
##               BoxID886               BoxID887               BoxID888 
##          -3.662041e-01          -3.662041e-01          -3.662041e-01 
##               BoxID889               BoxID890               BoxID891 
##          -3.662041e-01           7.202856e-01          -6.486367e-01 
##               BoxID893               BoxID894               BoxID895 
##          -8.796858e-01           1.564655e+00           8.253127e-01 
##               BoxID902               BoxID904               BoxID905 
##          -2.310491e-01          -2.310491e-01           2.310491e-01 
##               BoxID906               BoxID907               BoxID908 
##          -2.310491e-01           2.609308e-14                     NA 
##               BoxID910               BoxID911               BoxID912 
##           2.909693e-01           3.419921e+00           2.694443e+00 
##               BoxID913               BoxID914               BoxID915 
##          -5.966212e-01                     NA          -1.042869e+00 
##               BoxID916               BoxID917               BoxID918 
##           2.280182e+00           1.300658e+00           1.749008e+00 
##               BoxID919               BoxID920               BoxID921 
##           2.640270e+00           2.751076e+00           1.551320e+00 
##               BoxID922               BoxID923               BoxID924 
##           1.165503e+00           6.825643e-01           9.999452e-01 
##               BoxID925               BoxID926               BoxID927 
##                     NA           3.128042e-01          -1.735401e+00 
##               BoxID928               BoxID933               BoxID934 
##                     NA           1.589989e+00           7.336232e-01 
##               BoxID935               BoxID936               BoxID937 
##           1.042775e+00           1.260167e+00                     NA 
##               BoxID939               BoxID943               BoxID944 
##                     NA           1.516116e+00          -1.865386e-01 
##               BoxID945               BoxID946               BoxID947 
##          -3.662041e-01           1.021575e+00                     NA 
##               BoxID949               BoxID950               BoxID951 
##           5.660778e-01                     NA          -1.325746e-14 
##               BoxID952               BoxID953               BoxID954 
##           8.283022e-01          -7.636573e-14           6.931472e-01 
##               BoxID955               BoxID956               BoxID957 
##           6.931472e-01           7.415802e-15          -5.305915e-15 
##               BoxID958               BoxID959               BoxID960 
##           2.172013e-14                     NA                     NA 
##               BoxID962               BoxID963               BoxID965 
##           3.127565e-01                     NA           1.397662e-15 
##               BoxID966               BoxID967               BoxID968 
##           2.581443e-14           9.258647e-16           1.505527e-15 
##               BoxID969               BoxID970               BoxID971 
##           9.872099e-16           2.202295e-15           1.242385e-15 
##               BoxID972               BoxID973               BoxID974 
##                     NA           2.679866e-14           2.597101e-14 
##               BoxID975               BoxID976               BoxID977 
##           2.310491e-01           1.229626e+00           2.310491e-01 
##               BoxID978               BoxID979               BoxID980 
##                     NA           3.662041e-01           7.992984e-01 
##               BoxID983               BoxID984               BoxID985 
##          -6.931472e-01          -1.670937e-15                     NA 
##               BoxID986               BoxID987               BoxID988 
##           5.364793e-01           5.633805e-15           2.542049e-15 
##               BoxID989               BoxID990               BoxID991 
##                     NA           1.351550e-01           2.612920e-14 
##               BoxID992               BoxID993               BoxID994 
##                     NA                     NA          -1.425555e+00 
##               BoxID997               BoxID998               BoxID999 
##                     NA                     NA                     NA 
##              BoxID1001              BoxID1002              BoxID1003 
##           2.293092e+00           3.105316e+00           3.105194e-01 
##              BoxID1004              BoxID1005              BoxID1006 
##           4.620981e-01           1.675974e+00           2.020665e+00 
##              BoxID1008              BoxID1009              BoxID1010 
##           1.338101e+00           1.890225e+00                     NA 
##              BoxID1011              BoxID1012              BoxID1013 
##           1.843143e+00           2.305241e+00           4.783615e-01 
##              BoxID1014              BoxID1015              BoxID1016 
##           1.872661e+00          -8.581729e-01          -7.067545e-01 
##              BoxID1017              BoxID1018              BoxID1019 
##          -1.259497e+00          -1.586440e+00          -1.072959e+00 
##              BoxID1020              BoxID1021              BoxID1022 
##          -2.318848e+00           5.408944e-01                     NA 
##              BoxID1024              BoxID1025              BoxID1026 
##          -1.257647e-01           4.791165e-01          -1.856750e+00 
##              BoxID1027              BoxID1028              BoxID1029 
##          -1.625701e+00          -9.431370e-01          -6.168667e-01 
##              BoxID1030              BoxID1036              BoxID1037 
##                     NA           1.599639e+00          -2.902761e-01 
##              BoxID1038              BoxID1039              BoxID1040 
##           9.744082e-02           9.219137e-02           5.501486e-01 
##              BoxID1041              BoxID1042              BoxID1043 
##           8.073058e-01          -7.916353e-01           1.321561e+00 
##              BoxID1047              BoxID1048              BoxID1049 
##           1.093229e+00          -1.288854e+00           4.866554e-01 
##              BoxID1050              BoxID1051              BoxID1052 
##                     NA          -8.658211e-01          -1.544473e+00 
##              BoxID1056              BoxID1057              BoxID1058 
##           3.297124e-01           1.878840e-02                     NA 
##              BoxID1059              BoxID1060              BoxID1061 
##           2.259599e-01          -7.992984e-01          -2.402620e+00 
##              BoxID1062              BoxID1063              BoxID1064 
##           3.888723e-01          -3.264496e-01          -1.589096e+00 
##              BoxID1065              BoxID1066              BoxID1067 
##           8.988786e-02          -4.156094e-01          -1.336534e+00 
##              BoxID1068              BoxID1069              BoxID1070 
##          -3.261881e-01           3.135752e-02           2.849739e-01 
##              BoxID1071              BoxID1072              BoxID1073 
##                     NA          -9.217362e-01          -5.917434e-01 
##              BoxID1074              BoxID1075              BoxID1076 
##          -1.285542e-01          -5.813427e-01          -5.013591e-01 
##              BoxID1078              BoxID1079              BoxID1080 
##          -4.576830e-01          -2.097190e+00                     NA 
##              BoxID1081              BoxID1085              BoxID1086 
##          -2.839794e-01           6.433033e-01           2.720691e-01 
##              BoxID1087              BoxID1088              BoxID1089 
##          -3.945900e-01          -8.457582e-01           8.222468e-02 
##              BoxID1090              BoxID1091              BoxID1092 
##                     NA          -2.467177e-01          -5.527427e-01 
##              BoxID1093              BoxID1094              BoxID1095 
##                     NA                     NA          -8.257777e-01 
##              BoxID1096              BoxID1097              BoxID1098 
##          -7.405664e-01          -1.885996e+00          -7.402826e-02 
##              BoxID1101              BoxID1102              BoxID1103 
##          -2.399478e+00          -1.775544e+00          -2.992412e-01 
##              BoxID1104              BoxID1105              BoxID1107 
##          -4.734736e-02          -1.002926e+00          -3.658512e-01 
##              BoxID1108              BoxID1109              BoxID1110 
##           5.870310e-01          -1.371029e+00                     NA 
##              BoxID1112              BoxID1115              BoxID1116 
##           2.114275e+00          -2.299762e-02           1.199895e+00 
##              BoxID1117              BoxID1118              BoxID1119 
##          -1.777662e-01          -7.438118e-02           1.729313e-01 
##              BoxID1120            SA0:IndicG0            SA1:IndicG0 
##                     NA           2.660122e-01                     NA
vcov(m3)[420:423,420:423]
##                BoxID1119 BoxID1120  SA0:IndicG0 SA1:IndicG0
## BoxID1119   4.922688e-01        NA 5.750740e-17          NA
## BoxID1120             NA        NA           NA          NA
## SA0:IndicG0 5.750740e-17        NA 1.356993e-02          NA
## SA1:IndicG0           NA        NA           NA          NA
test <- vcov(m3)[418:423,418:423]
sqrt(test)
##                BoxID1117    BoxID1118    BoxID1119 BoxID1120  SA0:IndicG0
## BoxID1117   7.016187e-01 4.961194e-01 4.961194e-01        NA 7.476903e-09
## BoxID1118   4.961194e-01 7.016187e-01 4.961194e-01        NA 7.679951e-09
## BoxID1119   4.961194e-01 4.961194e-01 7.016187e-01        NA 7.583363e-09
## BoxID1120             NA           NA           NA        NA           NA
## SA0:IndicG0 7.476903e-09 7.679951e-09 7.583363e-09        NA 1.164900e-01
## SA1:IndicG0           NA           NA           NA        NA           NA
##             SA1:IndicG0
## BoxID1117            NA
## BoxID1118            NA
## BoxID1119            NA
## BoxID1120            NA
## SA0:IndicG0          NA
## SA1:IndicG0          NA
#se SA0IndicG0 = 0.11936296 avec lmer 
#se SA0:IndicG0 = 1.16490 avec lm (sans interaction)


#MYSELF
m3 <- lm(log(Nb_eggs+1) ~ pop_gen + hab_gen + SA*IndicG0 + SA  + BoxID, data = data_PREF_three)

#Estimates SA
cf <-coef(summary(m3,complete = TRUE)) 
indic <- grep("SA",rownames(cf))
SAcoef <- cf[indic,1]
#names(SAcoef) <- c("SAGen", "SANonGen")
SAcoef
##         SA1 SA1:IndicG0 
##   0.3842848  -0.2660122
#Estimates se
indic <- grep("SA",rownames(vcov(m3)))
seSAcoef <- sqrt(diag(vcov(m3)[indic,indic]))
#names(SAcoef) <- c("seSAGen", "seSANonGen")
seSAcoef
##         SA1 SA1:IndicG0 
##  0.07951066  0.11649003

5.2.3 Plot

(PLOT_pref_G0 <- plot_RTP_residuals(dataset = data_PREF_three, 
                                    trait = "Nb_eggs", gen = "G0"))

(PLOT_pref_G2 <- plot_RTP_residuals(dataset = data_PREF_three, 
                                    trait = "Nb_eggs", gen = "G2"))

(PLOT_GEN_pref_G0 <- plot_Genetic_Nongenetic_residuals(dataset = data_PREF_three, 
                                    trait = "Nb_eggs", effect = "Non-genetic"))

(PLOT_GEN_pref_G2 <- plot_Genetic_Nongenetic_residuals(dataset = data_PREF_three, 
                                    trait = "Nb_eggs", effect = "Genetic"))

6 PLOT LOCAL ADAPTATION

legend <- lemon::g_legend(PLOT_pref_G0)
# 

# 
# ## ALL GENERATIONS 
# LOCAL_ADAPTATION_PLOT <- cowplot::ggdraw() + 
#   cowplot::draw_plot(PLOT_pref_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.01, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_eggs_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.31, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_rate_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.61, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(legend, x = 0.93, y = 0.5, width = 0.0001, height = 0.0001) + 
#   cowplot::draw_plot(PLOT_pref_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.01, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_eggs_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.31, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_rate_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.61, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot_label(c("Generation: G0","Generation: G0","Generation: G1", "A", "B", "C", " ",
#                     "Generation: G2","Generation: G2","Generation: G3", "D", "E", "F", " "),  
#                   x = c(0.1,0.4,0.7, 0.01, 0.30, 0.61, 0.92, 0.10,0.4,0.7, 0.01, 0.3, 0.61, 0.92), 
#                   y = c(1,1,1, 0.98, 0.98, 0.98, 0.98, 0.5, 0.5, 0.5, 0.48, 0.48, 0.48, 0.48), 
#                   hjust = c(0,0,0,0,0,0,0,0,0,0,0,0,0,0), 
#                   vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
#                   size = 16) 
#  LOCAL_ADAPTATION_PLOT



## ALL GENERATIONS 
LOCAL_ADAPTATION_PLOT <- cowplot::ggdraw() + 
  cowplot::draw_plot(PLOT_pref_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_eggs_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_rate_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.01, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_pref_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_eggs_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_rate_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(legend, x = 0.915, y = 0.5, width = 0.000001, height = 0.000001) + 
  cowplot::draw_plot_label(c("Generation: G0/G1","Generation: G2/G3","A", "B", "C",
                    "D", "E", "F"),  
                  x = c(0.12,0.55, 0, 0.44, 0, 0.44, 0, 0.44 ), 
                  y = c(1,1, 0.99, 0.99, 0.66,  0.66, 0.33, 0.33), 
                  hjust = c(0,0,0,0,0,0,0,0), 
                  vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
                  size = 18) 
 LOCAL_ADAPTATION_PLOT

 cowplot::save_plot(file =here::here("figures", "Supplements_FigS4_Local_Adaptation_Residuals.pdf"),
                   LOCAL_ADAPTATION_PLOT,
                   base_height = 28/cm(1), base_width = 22/cm(1), dpi = 610)

7 PLOT GENETIC vs PLASTIC

legend <- lemon::g_legend(PLOT_pref_G0)

# Genetic_NonGenetic_PLOT <- cowplot::ggdraw() + 
#   cowplot::draw_plot(PLOT_GEN_pref_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.01, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_GEN_eggs_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.31, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_GEN_rate_G2 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.61, y = 0.5, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(legend, x = 0.93, y = 0.5, width = 0.0001, height = 0.0001) + 
#   cowplot::draw_plot(PLOT_GEN_pref_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.01, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_GEN_eggs_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.31, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot(PLOT_GEN_rate_G0 + theme(legend.position = "none", 
#                                         plot.title = element_blank()), 
#             x = 0.61, y = 0, width = 0.25, height = 0.45) + 
#   cowplot::draw_plot_label(c("Genetic effects", "A", "B", "C", " ",
#                     "Plastic effects", "D", "E", "F", " "),  
#                   x = c(0.4, 0.01, 0.30, 0.61, 0.92, 0.4, 0.01, 0.3, 0.61, 0.92), 
#                   y = c(1, 0.98, 0.98, 0.98, 0.98, 0.5, 0.48, 0.48, 0.48, 0.48), 
#                   hjust = c(0,0,0,0,0,0,0,0,0,0), 
#                   vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
#                   size = 16) 
#  Genetic_NonGenetic_PLOT
 
 
# 
# 
# cowplot::save_plot(file =here::here("figures", "GENETIC_NONGENETIC_First_Third.pdf"),
#                   Genetic_NonGenetic_PLOT,
#                   base_height = 18/cm(1), base_width = 34/cm(1), dpi = 610)
# 
# 
# 

 
## ALL GENERATIONS 
Genetic_NonGenetic_PLOT <- cowplot::ggdraw() + 
  cowplot::draw_plot(PLOT_GEN_pref_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_GEN_eggs_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_GEN_rate_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.01, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_GEN_pref_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_GEN_eggs_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PLOT_GEN_rate_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(legend, x = 0.915, y = 0.5, width = 0.000001, height = 0.000001) + 
  cowplot::draw_plot_label(c("Genetic effects","Plastic effects","A", "B", "C",
                    "D", "E", "F"),  
                  x = c(0.12,0.55, 0, 0.44, 0, 0.44, 0, 0.44 ), 
                  y = c(1,1, 0.99, 0.99, 0.66,  0.66, 0.33, 0.33), 
                  hjust = c(0,0,0,0,0,0,0,0), 
                  vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
                  size = 18) 
 Genetic_NonGenetic_PLOT

 #
 cowplot::save_plot(file =here::here("figures", "Fig5_GENETIC_NONGENETIC_First_Third.pdf"),
                   Genetic_NonGenetic_PLOT,
                   base_height = 28/cm(1), base_width = 22/cm(1), dpi = 610)

8 PLOT POP EFFECT - HETEROGENEITY

# library("dplyr")
PPMR_ALL_Pref_Stim_G0 <- plot_realdata(dataset = data_PERF, 
                                        trait="Nb_eggs" ,
                                        printcor = TRUE,
                                        gen = "G0" , 
                                        xlim = c(0,94), 
                                        ylim = c(0,94), 
                                        formula_Blanquart="log(Nb_eggs+1) ~  Test_environment + Population + SA + Test_environment:Original_environment",
                                        xaxis_labelprint = "Fecundity in original",
                                        yaxis_labelprint = "Fecundity in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
PPMR_ALL_Pref_Stim_G2 <- plot_realdata(dataset = data_PERF, 
                                        printcor = TRUE,
                                        trait="Nb_eggs" ,
                                        formula_Blanquart="log(Nb_eggs+1) ~  Test_environment + Population + SA + Test_environment:Original_environment",
                                        gen = "G2" , 
                                        xaxis_labelprint = "Fecundity in original",
                                        yaxis_labelprint = "Fecundity in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
## [1] "Populations that do not have measures in sympatry have been removed"
PPMR_ALL_Pref_G0 <- plot_realdata(dataset = data_PREF_three, 
                                        printcor = TRUE,
                                        trait="Nb_eggs" ,
                                        formula_Blanquart="log(Nb_eggs+1) ~ Test_environment + Population + SA + 
                    Test_environment:Original_environment + BoxID",
                                        gen = "G0" , 
                                        xaxis_labelprint = "Oviposition preference in original",
                                        yaxis_labelprint = "Oviposition preference in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
PPMR_ALL_Pref_G2 <- plot_realdata(dataset = data_PREF_three, 
                                        trait="Nb_eggs" ,
                                        formula_Blanquart="log(Nb_eggs+1) ~ Test_environment + Population + SA + 
                    Test_environment:Original_environment + BoxID",
                                        gen = "G2" ,
                                        printcor = TRUE,
                                        xlim = c(0,80), 
                                        ylim = c(0,80), 
                                        xaxis_labelprint = "Oviposition preference in original",
                                        yaxis_labelprint = "Oviposition preference in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
PPMR_ALL_Perf_G0 <- plot_realdata(dataset = data_PERF_Rate, 
                                         trait="Rate" ,
                                          gen = "G0" , 
                                        printcor = TRUE,
                                        xlim = c(0,1.1), 
                                        ylim = c(0,1.1), 
                                        formula_Blanquart="asin(sqrt(Rate)) ~  Test_environment + Population + SA + log(Nb_eggs) + Test_environment:Original_environment",
                                        xaxis_labelprint = "Offspring performance in original",
                                        yaxis_labelprint = "Offspring performance in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
PPMR_ALL_Perf_G2 <- plot_realdata(dataset = data_PERF_Rate, 
                                        trait="Rate" ,
                                        printcor = TRUE,
                                        gen = "G2" ,
                                  formula_Blanquart="asin(sqrt(Rate)) ~  Test_environment + Population + SA + log(Nb_eggs) + Test_environment:Original_environment",
                                        xlim = c(0,0.7), 
                                        ylim = c(0,0.7), 
                                        xaxis_labelprint = "Offspring performance in original",
                                        yaxis_labelprint = "Offspring performance in alternative")
## [1] "Converting test_environment column into a factor"
## [1] "Converting original_environment column into a factor"
## [1] "Populations that do not have measures in sympatry have been removed"
legend <- lemon::g_legend(PPMR_ALL_Perf_G2)

 
## ALL GENERATIONS 
PPMR_ALL <- cowplot::ggdraw() + 
  cowplot::draw_plot(PPMR_ALL_Pref_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PPMR_ALL_Pref_Stim_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.0, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PPMR_ALL_Perf_G0 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.01, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PPMR_ALL_Pref_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.66, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PPMR_ALL_Pref_Stim_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0.33, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(PPMR_ALL_Perf_G2 + theme(legend.position = "none", 
                                        plot.title = element_blank()), 
            x = 0.44, y = 0, width = 0.4, height = 0.3) + 
  cowplot::draw_plot(legend, x = 0.915, y = 0.5, width = 0.000001, height = 0.000001) + 
  cowplot::draw_plot_label(c("Generation G0/G1","Generation G2/G3","A", "B", "C",
                    "D", "E", "F"),  
                  x = c(0.12,0.55, 0, 0.44, 0, 0.44, 0, 0.44 ), 
                  y = c(1,1, 0.99, 0.99, 0.66,  0.66, 0.33, 0.33), 
                  hjust = c(0,0,0,0,0,0,0,0), 
                  vjust = c(1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5, 1.5),
                  size = 18) 
 PPMR_ALL

# 
# cowplot::save_plot(file =here::here("figures", "Fig4_POP_effect_ALL.pdf"),
#                   PPMR_ALL,
#                   base_height = 28/cm(1), base_width = 22/cm(1), dpi = 610)